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Abstract Analyses of two-particle correlations have provided the chief means for determining 
spatio-temporal characteristics of relativistic heavy ion collisions. We discuss the theoretical 
formalism behind these studies and the experimental methods used in carrying them out. Recent 
results from RHIC are put into context in a systematic review of correlation measurements 
performed over the past two decades. The current understanding of these results is discussed in 
terms of model comparisons and overall trends. 
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1 Introduction 

The study of nucleus-nucleus collisions at ultra-relativistic energies aims to char- 
acterize the dynamical processes by which matter at extreme densities is produced 
and the fundamental properties that this matter exhibits. In nucleus-nucleus col- 
lisions, how do partonic equilibration processes proceed? For how long, over 
which spatial extension, and at which density is a QCD equilibration state ap- 
proached, and what are its properties? Particle densities attained during a heavy 
ion collision are expected to exceed significantly the inverse volume of a hadron. 
This implies that the high temperature phase of QCD, the Quark Gluon Plasma, 
comes within experimental reach. Chiral symmetry restoration and deconfine- 
ment phase transition are testable in heavy ion collisions. However, the experi- 
mental study of QCD at high temperatures and densities is complicated by the 
short lifetime and mesoscopic extension of the produced system. Femtoscopy, the 
spatio-temporal characterization of the collision region on the femtometer scale, 
is needed to frame any discussion of dynamical equilibration processes. 

The Relativistic Heavy Ion Collider (RHIC) just completed the first part of a 
dedicated experimental heavy ion program. Center of mass energies (-y/sjvjv = 200 
GeV) exceeded those of previous fixed target experiments by a factor 10. The 
current discussion of RHIC data focuses mainly on several qualitatively novel 
phenomena that all support the picture that dense and rapidly equilibrating 
QCD matter is produced in the collision region (1; 2; 3; 4). In particular, identi- 
fied single inclusive hadron spectra appear to emerge from a common flow field 
whose size and dependence on transverse momentum and azimuth is consistent 
with expectations that the produced matter is a locally equilibrated, almost ideal 
fluid of very small viscosity. Moreover, high-pT hadron spectra show a centrality 
dependent, strong suppression in Au+Au collisions, but not in a d+Au control 
experiment, indicating that even the hardest partons produced in the collision 
participate significantly in equilibration processes. In short, experiments at RHIC 
have demonstrated already that heavy ion collisions produce dense and equilibrat- 
ing matter, and that controlled experimentation of this matter is possible using 
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a large variety of probes. 

Despite these successes, numerous questions remain concerning the state of the 
matter produced in these collisions. Most notably, the equation of state is far from 
being determined, and issues concerning chiral symmetry restoration are largely 
unresolved. Addressing these fundamental questions about bulk matter requires 
a detailed understanding of the dynamics and chemistry of the collision, which 
can only be acquired by thorough and coordinated analyses of data and theory. 
In particular, spatio-temporal aspects of the reaction need to be experimentally 
addressed. The small size, ~ 10~ 14 m, and transient nature, ~ 10~ 22 seconds, 
of the reactions preclude direct measurement of times or positions. Instead, 
femtoscopy must exploit measurements of asymptotic momenta. Correlations of 
two final-state particles at small relative momentum provide the most direct link 
to the size and lifetime of subatomic reactions (5; 6; 7; 8; 9; 10; 11; 12). Because 
correlations from either interactions or from identity interference are stronger for 
smaller separations in space-time, spatio-temporal information can be most easily 
extracted for small sources, unlike the limitations of microscopes and telescopes. 

The interference of two particles emitted from chaotic sources was first ap- 
plied by Hanbury-Brown and Twiss (13; 14), where photons were exploited to 
determine source sizes for both laboratory and stellar sources in the 1950s and 
1960s. Correlations of identical pions were shown to be sensitive to source di- 
mensions in proton-antiproton collisions by Goldhaber, Goldhaber, Lee, and Pais 
in 1960 (15). In the 1970s, these methods were refined by Kopylov and Pod- 
goretsky (16; 17; 18; 19), Koonin, (20), and Gyulassy (21), and other classes 
of correlations were shown to be useful for source-size measurements, such as 
strong and Coulomb interactions. Bevalac analyses showed that interferome- 
try was truly capable of quantitatively determining spatial and temporal source 
dimensions (22; 23) and providing a stringent test of dynamical models (24). 
Throughout the last 25 years this phenomenology has developed into a precision 
tool for heavy ion collisions. Whereas hadronic sources are short lived and one 
measures correlations of the momentum of outgoing particles, stars are long lived 
and require experimental filters to enforce the approximate simultaneity of the 
two photons. Although the theory for these two classes of measurement are very 
different (25; 26), the heavy-ion community often uses the term HBT, in refer- 
ence to Hanbury-Brown and Twiss's original work with photons, to refer to any 
type of analysis related to size and shape that uses two particles at small relative 
momentum. To some, however, the term HBT refers only to identical-particle 
interferometry. Following Lednicky, we will employ the term femtoscopy (27; 28) 
to denote any measurement that provides spatio-temporal information, including 
coalescence analyses. 

Femtoscopic measurements from truly relativistic heavy ion collisions were first 
reported almost twenty years ago. Since then, measurements have been per- 
formed for collisions at energies spanning two orders of magnitude. In a double 
sense, then, this review examines recent RHIC results within the larger context of 
two decades' worth of femtoscopy. The theory and phenomenology of correlation 
femtoscopy are reviewed in the next section, with particular emphasis on describ- 
ing the approximations used to derive the connection between spatio-temporal 
aspects of the emission function and correlations constructed from final-state mo- 
menta. Experimental methods and techniques are correspondingly reviewed in 
Section 3. Section 4 presents experimental results, with an emphasis on describing 
the systematics of source dimensions and lifetimes as a function of beam energy, 
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system size, particle species and a particle's momentum. In addition to source 
dimensions, results for phase space density and entropy are presented. Compar- 
isons of experimental results and transport models are presented in Section 5, 
with an emphasis on explaining the "HBT puzzle," i.e., the fact that dynamic 
descriptions that incorporate a phase transition to a new state of matter with 
many degrees of freedom significantly over-predict observed source sizes. Section 
6 summarizes the current state of the field and lists new directions and challenges 
for future theoretical and experimental analyses. 



2 Theory and Phenomenology Basics 
2.1 Formalism 

Two-particle correlation functions are constructed as the ratio of the measured 
two-particle inclusive and single-particle inclusive spectra, 

rab(J} , dN ab /(d 5 p a d 3 p b ) 

{ ' qj {dN a /d?p a )(dN b /d 3 p b )' [) 

P - Pa +Pb , ^=<g^- fr'^- i V 

The theoretical analysis of (1) aims at relating this experimentally measured 
correlation to the space-time structure of the particle emitting source (7; 6; 8; 11). 
Two forms are common for connecting the measured correlation function to the 
space-time emission function s(p, x) through a convolution with the wave function 
<f>. In the first form (29), 

C ab, p x = I d 4 x a d 4 x b s a (p a , x a )s b (p b , x b )\<j>(<j, r')| 2 
J d 4 x a s a (p a ,x a ) J d 4 x b s b (p b ,x b ) 

In calculations of two-particle correlation functions, the squared relative two- 
particle wave function \(f)\ 2 serves generally as a weight, and the emission function 
s(p,x) contains all space-time information about the source because it describes 
the probability of emitting a particle with momentum p from a space-time point 
x. Here, and throughout this section, primes denote quantities in the center-of- 
mass frame, i.e., the frame where P = 0. The source function s a is evaluated at 
the momentum p a = m a P/ (m a + m b ), p° a = E a (p a ). 

The second form, which is equally justified as Equation 2 by the approximations 
described further below, is, 

C a6 (P,q)-l = |dV<S P (r')[|<Kq'Y)| 2 -l], (3) 
J d 4 x a d 4 x b s a (p a , x a )s b (p b , x b )5 (V - x' a + x! b ) 



Sp(r') = 



f d 4 X a d 4 X b S a (p a , X a )s b (p b , X b ) 



This expression allows one to consider \cj)\ 2 as a kernel with which one can trans- 
form from the coordinate-space basis to the relative-momentum basis. It also 
emphasizes the limitation of correlation functions, that they can provide, at best, 
the function <Sp(r'), the distribution of relative positions of particles with iden- 
tical velocities and total momentum P as they move in their asymptotic state. 
Thus, correlations do not measure the size of the entire source. Instead, they ad- 
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dress the dimensions of the "region of homogeneity," a term coined by Sinyukov 
(30), i.e., the size and shape of the phase space cloud of outgoing particles whose 
velocities have a specific magnitude and direction. If the collective expansion 
of the produced matter is strong, as is the case in central collisions, then the 
region of homogeneity is significantly smaller than the entire source volume. In 
the following, we discuss the assumptions on which Equations 2 and 3 are based. 

We start from explicit expressions for the one- and two-particle spectra in terms 
of T-matrix elements. For one-particle emission, 



E^- = J d 4 x s(p,x) = YL\j d 4 xT F ,(x)e- ip - x \ 2 , 



(4) 



a(p, x ) = Yl I d A 5xT* F ,{x + 5x/2)T F ,{x-5x/2)e iv - Sx . (5) 

Here, F' refers to the state of all other particles in the system. All interactions 
with the residual system are incorporated into the T matrix. However, there is a 
choice as to whether mean-field interactions are included in the T matrix or are 
instead incorporated into the evolution matrix (31; 32; 33; 34; 35). For example, 
one can include the Coulomb interactions with the residual system by replacing 
the phase factor e ip ' x in Equation 4 with an outgoing Coulomb wave function. 
This can be quantitatively important, in particular for slow particles. It becomes 
difficult when the two particles interact with one another through the strong or 
Coulomb force, as this represents a quantum three-body problem. 

Assumption 1: Higher order (anti)symmetrization can be neglected. Equa- 
tion 4 implies that all particles with asymptotic momentum p must have had 
their last interaction with the source at some point x. For distinguishable par- 
ticles, this is indeed the case and Equation 4 does not represent an assumption. 
However, if there are N a > 1 particles of the same type a, then one must consider 
T a (x\ ■ • • x jv a ) ■ The evolution matrix is then no longer a simple phase factor but 
includes N a \ interference terms. The single-particle probability is then obtained 
by integrating over the other iV a — 1 momenta. This can be accomplished ex- 
plicitly for simple source functions (36; 37; 38; 39; 40). The distortion to the 
single-particle spectra and to the two-particle correlation function were found 
to be important when the phase space density approached unity. Otherwise, 
Equation 4 is well justified. 

Assumption 2: The emission process is initially uncorrelated. In writing Equa- 
tion 4, one requires that two-particle matrix elements factorize, Tpn(x a , Xf,) = 
Tp^x^Tp^Xb), i.e., that the emission is independent. If multi-particle sym- 
metrization can be neglected, the two-particle evolution operator factorizes into 
a center-of-mass and a relative operator. One has U(x a , Xb',p a ,Pb) = u q'( x ' a ~ 
x' b )e lP '( EaXa / Minv+E b Xb / Mmv } for non-identical particles, whereas for identical par- 
ticles U = e lP ^ Xa+Xb ^ 2 {u q i{x' a - x' b ) ± {u q i(x' b - x' a )) / \ / 2. This is illustrated in 
Figure 1. Then, the two-particle probability can be expressed in terms of one- 
particle source functions, 



EaEh Jp&p h = J dix ^bd 4 qs a {{E'JM inY )P + q, x a )s b {(E' h /M inv )P - q, 



Xb) 



xd 4 dr'e^ 6r 'u* ql (x' a -x' b + 5r' /2)u q ,(x' a - x' b - Sr'/2). (6) 
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if a identical to b 




Figure 1: Schematic representation of the squared emission amplitude for two 
particles emitted independently from the grey-shaded source region and interact- 
ing with each other in the final states. For identical bosons (+) and fermions (-), 
correlation also involves interference between the paths. 



Assumption 3: Smoothness approximation (41; 42; 43). Equation 6 is difficult 
to evaluate as it requires knowledge of the source function evaluated off-shell. For 
the special case where the particles do not interact aside from identical particle 
interference, u q (x a — x b ) = [e tq '( Xa ~ Xb > ± e iq '( Xb ~ Xa ']/\/2, the integrals over q and 
Sr' can be performed analytically, 



The source functions in the interference term are evaluated off-shell for non-zero 
q, Pq/2 ^ E(P/2). The smoothness approximation replaces s(P/2, x a )s(P/2, x b ) 
with either s{E(P/2),P/2,x a )s(E(P/2),P/2,x b ), which leads to Equation 3, or 
with s(p a , x a )s(p b , x b ), which leads to Equation 2. If the first approximation is 
performed, one should also make the same approximation for the denominator. 
The smoothness approximation has been checked for expanding thermal sources, 
where it was found to be very reasonable for large (RHIC-like) sources, but quite 
questionable for smaller sources such as those found in pp or e + e~ collisions (43). 

Assumption 4 : Equal time approximation. For the general case where the evo- 
lution operator incorporates Coulomb or strong interactions, deriving Equations 3 
and 2 from Equation 6 is more complicated. First, the smoothness assumption 
amounts here to neglecting the q dependence in the product of the source func- 
tions in Equation 6. This assumption is somewhat more stringent in the presence 
of final state interactions because the relevant range of q extends beyond q. With 
this assumption, one obtains a (5-function constraint for 5r' , and the integrand 
of (6) is proportional to the squared evolution matrix \u' q (x' a — x' b )\ 2 . This evo- 
lution matrix has non-zero time components, which must be neglected if one is 
to identify it with the relative wave function, because the relative motion in the 
pair rest frame is small, one expects this approximation to be reasonable, but it 
has not yet been tested in model studies. 

The above formalism is semi-classical in the sense that a quantum-mechanical 
particle emission probability, defined by the T-matrix elements, is usually ap- 
proximated by classical source functions. As a consequence, quantum uncertainty 





±s(P/2, X a )s(P/2, X b ) COs[(> a - p b ) ■ (x a ~ X b )]} . 
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limits the applicability of Equations 3 and 2. To illustrate this limitation, source 
functions have been evaluated by convoluting the emission function with wave 
packets(44; 45; 46; 47) of spatial width a. This leads to a broadening of spatial 
distributions by (AR) 2 ~ 0.5 — 1.0 fm 2 . because quantum smearing may already 
be incorporated into some of the semi-classical treatments through the choice of 
the initial density distribution, these calculations should be regarded as indicative 
of the theoretical uncertainty, because this error affects the size in quadrature, 
it is negligible for large sources, but might be significant for sources near 1 fm 
in size with strong space-time correlations (42; 43; 45; 46; 48). In particular, 
analyses of tttt correlations from e + e~ collisions, which result in source sizes of 
less than 1.0 fm(49; 50), are difficult to interpret in the above formalism. 



2.2 Identical-Particle Interference 

In the absence of strong and electromagnetic final state interactions, the wave 
function of an identical particle pair in Equation 3 becomes 

|0(q',r')| 2 - 1 = ±cos(2q'-r'). (8) 

The distribution of separations in coordinate space <Sp(r') can then be determined 
by performing a three-dimensional Fourier transform of C(q') — 1 . It is instructive 
to consider the properties of this inversion in more detail. The curvature of C(q) 
at q = can be related to the mean-square separation of the three-dimensional 
shape of Sp(r) (we neglect the P labels in C and S), 

= (nrj) = [ d 3 rS(r)r irj . (9) 

9=0 J 

This relation has been useful to qualitatively illustrate the relation between spe- 
cific space-time information and specific features of the correlator. However, 
applying the identity quantitatively requires careful consideration of pions from 
longer-lived resonances which can dominate the calculation of (r 2 ) if not ac- 
counted for. 



2.3 Correlations from Coulomb and Strong Interactions 

Compared to the case for non-interacting identical particles, where the transfor- 
mation between <Sp(r) and C(P, q) is a Fourier transform, analyzing the experi- 
mentally measured correlation function with Equation 3 to determine the source 
function is more complicated. Understanding the resolving power of the kernel 
\4>(q_, r| 2 requires a detailed understanding of the relative wave function. Once one 
averages over spins, the squared relative wave function is a function of q, r and 
cos6 qr . In relativistic collisions, correlation analyses are usually confined to light, 
singly charged hadrons. Coulomb-induced correlations are then weak and must 
be analyzed at small q, where quantum effects become important (qr/h ~ 1). 
The relative two-particle wave function in the presence of Coulomb interactions 
can then be written as a function of qr/h, r/ao and cos0 gr . 

4> = r(l + ir ? ) e -^/ 2 e^ r |l + f;/ l „.(r/aor|, (10) 
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where ao is the Bohr radius, hi = 1 and h n = n _] n ^ v h n -\. Here, rj = fie 2 /hq 
is independent of r, and for small r/ao the correlation function behaves as the 
Gamow factor, G{rf) = e _7 "»|r(l + ir])\ 2 = 2tttj / (e 27rr > - 1). Thus, the Coulomb 
kernels have little resolving power for nir correlations where ao = 387 fm, but 
have greater resolving power for pK correlations where ao = 84 fm. 

For r << ao, the effects of Coulomb interactions can be removed easily from 
the correlation function because r\ is independent of r (21), 

|0| 2 « G(rj)[l + cos(2grcos0 ?r )] . (11) 

For realistic source sizes, the order r/ao corrections are of the order of 10% 
for 7T7T correlations and are larger for heavier pairs (41; 51). Significant effort 
has been invested in "correcting" experimental correlation functions to remove 
Coulomb effects to all orders, but such corrections are model-dependent. The 
safest method for determining the source function is either inverting the full kernel 
(52; 53; 54; 55; 56; 57) or fitting C(q) to some parameterized form for S(r), which 
is convoluted with the full kernel. Neither of these tactics are computationally 
prohibitive. 

Strong interactions can also be exploited to provide size and shape information. 
If the size of the source is much larger than the range of the potential between 
the two particles, the kernel (\<p\ 2 — 1) can be determined entirely from knowledge 
of the phase shifts. Pairs that have a resonant interaction are especially useful, 
because the resonance will lead to a peak whose height is inversely proportional 
to the source volume, if qR » h. At small q the kernel is determined by the 
scattering length, a (58), and the height of the correlation at q = becomes 

C(q = 0)-l = (-^ + £), (12) 

where the averaging is performed using <Sp(r) as a weight. Of course, the effects 
of strong interactions, Coulomb interference, and identical-particle interference 
can all combine as is the case for pp correlations. The pp kernel has been analyzed 
in depth by Brown and Danielewicz, where the kernel was inverted and applied to 
experimental pp data. Evidence was seen for significant non-Gaussian behavior in 
the resulting source functions (56; 57). Strong interactions also provide angular 
resolving power (59) which can be understood from the perspective of classical 
trajectories. Even s-wave scattering can be exploited to discern information about 
shape. 

Strong and Coulomb-induced correlations apply to both identical and non- 
identical particles. For non-identical particles, the wave function is not sym- 
metrized and \<p(n, r)\ 2 ^ |</>(q, — r)| 2 , which results in odd components of the 
correlations function, C(q) / C(— q) if there are odd components of S(r) as 
is the case for non-identical particle pairs. This asymmetry can be experimen- 
tally exploited to investigate the spatio-temporal differences between the emission 
functions of different particle species (60; 61); this requires, however, sufficient 
statistics to select on the angle between the total and relative momentum in the 
pair center of mass Z(q, P) (62). 
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2.4 Coordinate Systems 

Correlation functions depend on two three-dimensional momenta, P and q. For 
high-energy collisions, analyses are usually performed in the longitudinally co- 
moving system (LCMS), a rest frame moving along the longitudinal (beam) di- 
rection such that P z = 0. Axes are usually chosen according to the out-side-long 
prescription. The longitudinal axis is chosen parallel to the beam, while the 
outward axis points in the direction of P, which is perpendicular to the beam 
axis. The sideward axis is chosen perpendicular to the other two. If the system 
is boost-invariant, observables expressed in the LCMS variables are independent 
of P z and the source has a reflection symmetry about the ri ong = plane. If 
the collision is central, there is also a reflection symmetry about the r s jd e = 
plane. Any four-vector V can be expressed in this coordinate system using the 
four-momentum P to project out the components(63; 64; 65), 

Mong = (P V Z -P Z V )/M T , (13) 

V Qut = (P X V X + PyVy)/P T , 
V side = (P X Vy-PyV X )/P T , 

where Mj, = Pq — P 2 and P^, = P 2 + P 2 . Dimensions of the source function are 
typically quoted in this coordinate system. One could also perform a second boost 
to the pair frame, in which the transverse components of the total momentum 
are zero. Then, 

, _ M hrv (P x V x +PyVy) P T , , 

out M T P T M T M inv ' { } 

where M? nv = P 2 . Relative wave functions are more conveniently expressed in 
the frame of the pair. For instance, a sharp resonant peak is no longer sharp if the 
correlation is viewed away from the pair frame. For pairs where the correlation 
is influenced by Coulomb and strong interactions, most analyses are conducted 
in the pair frame. 

For non-zero impact parameters, azimuthal symmetry is lost and source func- 
tions also depend on the azimuthal direction of the pair's total momentum (66; 
67; 68). Also, if boost-invariance is broken, the pair's rapidity needs to be speci- 
fied. In this more general case, reflection symmetries are broken and the choice of 
the coordinate axes becomes somewhat arbitrary. One could orient the 
cording the event's impact parameter, or one could rotate the coordinate system 
so that in the new frame cross terms such as (xy) vanish (illustrated in Figure 2). 
One would then specify the Euler angles as part of the description of the source 
function. 



2.5 Gaussian Parameterizations 

To gain a physical understanding of the three-dimensional spatio-temporal source 
distribution, it is useful to summarize its size and shape with a few parameters. 
This motivates the study of Gaussian parameterizations for the source <Sp(r') 
and the two-particle correlator. Realistic sources deviate from Gaussians, e.g. 
by exponential tails caused by resonance decay contributions. The extracted 
Gaussian source parameters may thus depend on details of the fitting procedures. 
These shortcomings can be overcome with imaging methods, or more complicated 
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forms for the fitting. A more general three-dimensional analysis of correlations 
would involve decomposing both the correlation and source functions in either 
spherical or Cartesian harmonics (69; 70). Although the detailed non-Gaussian 
aspects of the correlation are important, the extra information can also cloud 
the main trends in the data. In practice, Gaussian parameterizations provide the 
standard minimal description of experimental data. 

2.5.1 The general case The most general form for a Gaussian source is 
exp{— A a fj{x a — x a )(xp— xp)}, where x a refers to the four dimensions x, y, z, t, and 
A is a 4 x 4 real symmetric matrix. The most general form has 14 parameters, 
10 parameters for A and four more parameters for the offsets x a . Reflection 
symmetries can be used to eliminate certain cross terms and some of the offsets 
(71). Furthermore, if the particles are identical or have the same phase space 
distributions, all the offsets can be set to zero, because the source function 
for the second species b might also have 14 parameters, there could be up to 
28 Gaussian parameters in describing both s a and Sb- However, because <S(r) 
depends only on the distribution of relative spatial coordinates in the pair frame, 
only nine Gaussian parameters are required to describe the most general S for a 
given P. Three of these nine parameters can be identified with Gaussian widths, 
three can be identified with offsets, and the last three can either be identified 
with cross terms or with the three Euler angles describing the orientation of the 
three-dimensional ellipse. 

Using the reflection symmetries for mid-rapidity sources in a symmetric central 
collision, a Gaussian parameterization of the emission function for particle species 
a in the out-side-long coordinate system reads 

/ \ I (^out ^a,out — V s a (t — ia)) /-, r \ 

s a (p,x) ~ exp<^ '—— 2 '■ (15) 

^ ZU a,out 

(^sidc) (-^long) (t ^a) 1 



2< si dc 2i% long 2(Ar a )2 j ■ 

The symmetries forbid any cross terms in the exponential involving x s ^ c or 
x\ ong , such as x S ide^out or x\ ong t, but do not forbid a cross term between outward 
position x out and time t. Here, this cross term is taken into account by allowing 
the source to move in the outward direction with a velocity V s . The symmetries 
also forbid offsets in the sideward or longitudinal direction. These other offsets 
have also been addressed for non-central collisions (72; 73). 

The correlation function is determined by the phase space density of the final 
state, Equation 27. The resulting phase space density is 



/o(p,r,t) ~ exp 



[x out - X a (t)] 2 xj ide zLg 



Vout + (Vs,a ~ V ± nAr a )2] 2R a sidc2 2R a long z j ' 
X a (t) = x a , ont + V ± (t-t a ). (16) 



Here, V± is the velocity of the pair in the LCMS frame. 

The correlation function is determined by the relative distance function iSp(r'), 
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see Equation 3, 

Sv(r') ~ cxd f ^° ut - ^° ut ^ ^ks ri ° n s \ (u) 

P{) "I 4 7 i^ ut 4^ ide 4^ ong j (17) 

#out = l[Rl out +Rl out + (V s , a -V ± ) 2 (AT a ) 2 + (V s , b -V ± )\AT b ) 2 ], 

E>2 _ I" r>2 _i_ r>2 1 c>2 _ _ I" r>2 , p2 1 

-'•'side r> L «,side ' JX 6,sidcJ > -'''long r> L a,!""" 



Jong ~r -"-bjongj > 



-^out — x a,out x b,out) 



where 7_i_ = (1 — Vj_) 1 ^ 2 - Thus, there are four measurable parameters, -Rout, 
-Rsido -Rlong) and ^out- in the absence of any symmetry, there are five more 
terms: three cross terms (^ ut;Side , i?o UtJong and #side ,iong)' an d two more offsets 
(X si( j e and -X"i ong ). For identical particles, all the offsets are zero. 

2.5.2 Sensitivity to Lifetime Given the symmetries used to derive Equa- 
tion 17, the most experiment can provide is the determination of the four param- 
eters, -Rout) -Rsido -Riong> and ^out- The only way to independently determine 
the lifetime is to assume that the two transverse spatial sizes are approximately 
equal (74). After assuming R 2 a out + R 2 out « R 2 sidc + i?f )Bide , Equation 17 yields 
for identical particles (V a , s = H,s = V s , Ar a = Ar b = At), 

(V ± -V s ) 2 (Ar) 2 ^R 2 out -R 2 sidc . (18) 

In general, however, the source velocity is not precisely known and outward and 
sideward spatial dimensions are not exactly equal; these factors result in a signif- 
icant systematic error when applying Equation 18, especially because temporal 
effects enter in quadrature. The -R ou t/-Rside ratio only provides a reliable estimate 
of the lifetime when (V± — V s )At is much larger than the transverse size. 

2.5.3 Gaussian Cross Terms In addition to the three Gaussian parame- 
ters, -Rout) -Rside and -Riong) that describe the dimensions of a Gaussian source, one 
needs, in the general non-symmetric case, three more parameters to describe the 
angular orientation of the principal axes. Figure 2 displays how the principal axes 
might differ from the out-side- long axes once the collisions are off center. These 
three Euler angles combined with three sizes can also be related to the parameters 
Aij in the general form for a three-dimensional Gaussian, exp(— Aij^rj), where 
A is a symmetric matrix with six independent components (75; 71; 66; 68). 

For pairs of identical particles, the correlation function for Gaussian sources is 
also Gaussian, 

C(q) = 1 + exp(-(^ 2 ) ij Q l Q,) , (19) 

where Q = 2q. The six experimentally determined parameters, (TZ) 2 j, can be 
related to the moments of S(r') (76), 

{nr 3 ) = (TZ 2 ) ir (20) 

For central collisions, the source sizes (R 2 )ij depend only on the longitudinal 
pair momentum Pl and on the modulus of the transverse pair momentum \Pt\- 
This is different for non-central collisions for which the azimuthal direction </> pa ir = 
Z(Pt, b) of the transverse pair momentum with respect to the impact parameter 
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^out-long 




Figure 2: For non-central collisions the principal axes describing the orientation 
of the region of homogeneity can differ from the out-side-long axes. By viewing 
the source distribution from the perspective where the beam axis is oriented 
horizontally (left panel) and from peering down the beam pipe (right panel), the 
orientations leading to out-long and out-side cross terms are illustrated. 

direction b matters. The </> pa i r -dependences of (R 2 ) Q ut, (R 2 )side, and (-R 2 ) ou t-sidc 
then characterize the degree to which the initially out-of-plane-extended source 
geometry has expanded to the point where it becomes in-plane-extended (66; 67). 
The out-longitudinal and side-longitudinal cross terms contain information about 
the extent to which the main axis of the emission ellipsoid is tilted with respect 
to the beam axis(67). We return to this topic in Section 4.2. 

The Yano-Koonin parameterization (77; 78; 71; 79) provides an alternative 
basis for describing the out-long cross term. The Yano-Koonin form is based on 
the assumption that one can boost along the beam axis to a source frame where 
the correlation function has a simple form. 

C(P, Q) = 1 + exp [-Q\r\ - Q\R\ - Q 2 R 2 ] , (21) 

where Q is the momentum difference defined in the source frame which has ra- 
pidity j/yk- In that frame Rq is the Gaussian lifetime and R± and Ri\ are the 
dimensions of the source. This can be transposed to the out-long-side frame by 
boosting Q along the beam axis to the frame where P z = 0, i.e. the LCMS 
frame, then using the fact that Qo = QoutP±/Po = QoutV± in the new frame. 
This yields Q = cosh(y 7r7r - yyK)QoutV± - sinh(y 7r7r - yviOQiong, and Q\\ = 
cosh(y 7r7r - yYK)Q\on g ~ sinh(y 7r7r - yyK)QoutV±, where y W7r is the rapidity of the 
LCMS frame. Substituting these expressions into Equation 21 yields a cross term 
in the exponential equal to QoutQiong(-Ro + smn ^{v^ ~ VYk), which disap- 
pears when yyx = ynw By fitting yyx as a function of the pair rapidity, aspects 
of boost invariance can be tested. Given that the distribution of source rapidities 
should fall off for large rapidities, one expects yyx to lag y W7r , because pions of 
a given rapidity would more likely have been emitted from sources with smaller 
rapidities (78). 

2.6 The A factor 

Many pions measured in experiment come from long-lived decays. Pions from 
weak decays may or may not, depending on the experiment, be identified and re- 
moved from the analysis, because their decay vertices are typically a few centime- 
ters from the reaction center. Decays from r] or n' mesons occur a few thousand 
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fm away from the center of the collision. At these distances, they are effectively 
uncorrected with other particles but cannot be identified with experiment. If 
a fraction A of the pairs originate from the spatio-temporal region relevant for 
correlations, the correlation is muted by the factor A (21). If the source function 
is divided into two contributions, <S(r) = A5i OC ai + (1 — A)5oo, where both <Si OC ai 
and Soo integrate to unity, the resulting correlation is 



If the experimental sample includes a contamination from weak decays, rj, or 
mis-identified particles of a fraction /, the lambda factor becomes (1 — f) 2 . It is 
not uncommon for this contamination factor to be near 30%, which results in A 
near one half. 

Certainly, this division is somewhat artificial, as there are non-Gaussian tails, 
or halos (80), to 5i OC ai due to such causes as the exponential fall-off of the source 
function in the longitudinal direction or semi-long-lived resonances such as the (ft, 
whose lifetime is 40 fm/c. Non-Gaussian behavior is a subject for either imaging 
(52; 53; 54; 55; 56; 57) or for more complicated parameterizations. 

The A factor is often referred to as an incoherence factor, the name being mo- 
tivated by the properties of a coherent state, exp{i f dpn(p)[a(p) + a'(p)]}|0), 
which for identical particles leads to no correlation. Coherent states represent 
highly correlated emissions caused by phase coherence and thus violate the as- 
sumption of incoherence or uncorrelated emission implied by Assumption 2 as 
described early in this section. The question of whether an observation of A / 1 
is due to coherence or due to contamination from particles from far outside the 
source volume can be tested by analyzing three-particle correlations (81; 82). 
Such analyses of data at both SPS and RHIC have been consistent with the in- 
coherence conjecture (83; 84; 85). Microscopic model calculations at the AGS 
reproduce the excitation function of A when resonances contributions are in- 
cluded (86). 

2.7 Collective Flow and Blast- Wave Models 

Both longitudinal and radial collective expansion reduce the size of the region of 
homogeneity, i.e., the relevant volume for particles of a given velocity. For an 
infinite volume, the size of this region is set by the length one must move before 
collective velocity overcomes the thermal velocity, R ~ Vtherm /(dv/dz). 
2.7.1 Longitudinal Flow The Gaussian lifetime At described in Equa- 
tion 18 represents the spread of the emission times. A small value of Ar does 
not imply that particles were emitted early, but that they were emitted suddenly. 
Inferring the mean time at which particles are emitted requires a different as- 
sumption. For instance, at RHIC, the initial nuclei are Lorentz contracted by 
a factor of 100, and if there were no subsequent expansion, -Ri ong would be less 
than a Fermi. If one assumes that the system expands along the beam axis with 
no longitudinal acceleration, the collective velocity becomes 



If emission then comes from sources moving over a large range of rapidities (a 
boost-invariant expansion), the dimension along the beam axis for the source 




(22) 



Vco\\,z = z/t. 



(23) 
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Figure 3: because particles with heavier masses have smaller thermal velocities, 
their source volumes are more strongly confined by collective flow. For longitudi- 
nal flow (left panel) this results in smaller values of i2i on g f° r particles with higher 

ttit = \J m 2 + p^,. For radial flow (right panel) this confines heavier particles 
toward the surface, which results in both a reduced volume and an offset Ar in 
the outward direction. 



emitting zero-rapidity particles is determined by the distance one can move be- 
fore the collective velocity overwhelms the thermal velocity to force the emission 
function back to zero. The size can then be expressed as: 



-^long ~ , / , — M,herin(^/ • 

dv/dz 



(24) 



Whereas Rout/Rside gives information about the suddenness of emission, i?i ong 
provides insight into the mean time at which emission occurs given an estimate 
of the thermal velocity. 

For a thermal source with relativistic motion, the thermal velocity along the 
beam axis is determined by the temperature and the transverse mass, tut = 

J m 2 + p\ (63). For large mj- the thermal velocity in the longitudinal direction 

becomes non-relativistic, Vtherm = y/T/rriT, and the source size falls as 1/^mr 
which is referred to as tut scaling (87). This is illustrated in Figure 3. However, 
this assumes all particles are emitted with the same Bjorken time tb and tem- 
perature, independent of the transverse mass, because particles with high wit 
are probably emitted at lower tb, and because the temperature roughly behaves 
at T B 4//3 , the longitudinal size could fall even more quickly than m T 1 ^ 2 . 

In a boost invariant expansion, emission is a function of the Bjorken time 

tb = Vt 2 — z 2 , not the time t, and because t = \Jt b + z 2 , those particles emitted 
with small z have a head start. This is sometimes referred to as an inside- 
outside cascade. The transverse shape of <S(r) is then affected non-trivially by the 
expansion along the beam. The resulting correlation function can be calculated 
analytically in the case of pure identical-particle correlations (88; 89). 

Boost invariance is incorporated into blast-wave models with transverse expan- 
sion and assumed for many hydrodynamic models. The finite size of the system 
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would alter the results for two reasons. First, if the distribution of sources cov- 
ers only a finite range in 77, the tails of the distribution S(r) are chopped off. 
Assuming the distribution in 77 is Gaussian rather than uniform, 

111 

#2 7/2 T T + n T T T ' (25) 

^long v therm 7 B 'lG T B 

where t)q is the range of rapidities over which the sources are distributed. If rjc 
were 1.5 units of rapidity, the extracted values of tb from boost-invariant pictures 
would be underestimating tb by ~ 10%. 

A second shortcoming of boost-invariant models is that they ignore acceleration 
in the longitudinal direction. Accounting for this acceleration would alter the 
relation between the time and the velocity gradient. Neglecting this acceleration 
could also lead to a modest underestimate of tb (180). 

2.7.2 Transverse Collective Flow because transverse collective flow 
is intimately related to the pressure and viscosity, it is of central importance. 
Blast-wave models are based on pictures of thermal sources superimposed onto 
the transverse and longitudinal collective velocity profiles. Simple forms are then 
chosen for the profiles. Only two parameters are important for analyzing spectra, 
the temperature and the transverse velocity, because heavier particles are more 
sensitive to flow than are light particles, the two parameters can be adjusted to 
fit the spectra of several species. 

In addition to the temperature and transverse velocity, correlation measure- 
ments are also sensitive to the space-time parameters of the blast wave. In a 
minimal parameterization this would include the lifetime tb and the transverse 
size R. More sophisticated models would also include a spread in lifetime At 
and a surface diffuseness AR. Additional parameters ensue when one considers 
sensitivity to the reaction plane. Then, two parameters are needed to describe 
the transverse size, and two parameters are required to describe the transverse 
collective velocity. These parameters can then depend on the azimuthal direction 
of P (90). 

Choosing a blast- wave parameterization involves a number of choices about the 
form of the parameterization (91). Chemical potentials and temperatures might 
be chosen to vary with the transverse position r (92; 93) or might be chosen 
to be uniform. A wide variety of parameterizations have been employed for the 
transverse velocity profile, which might choose linear profiles for either v, jv or 
the transverse rapidity asinh(^v). In some parameterizations, the velocity profile 
has been chosen to rise quadratically with r (94; 95). Although hydrodynamics 
has been invoked as justification for different parameterizations, profiles from 
hydrodynamics vary according to the equation of state. 

For particles moving much faster than the surface velocity, transverse flow 
manifests itself by constraining particles to an increasingly small fraction of the 
blast-wave volume for the same reason that i?i on g fells with mx owing to longi- 
tudinal expansion (96). For large my, this leads to both i? s id e and R out falling as 
VV™T (71; 92; 93; 97). The fact that transverse dimensions fall with my might 
also result from the dynamics of cooling, superimposed with a growing fireball. 
This correlates high-energy particles with earlier times when the fireball was both 
smaller and at a higher temperature. 

Non-identical particles are of special interest in a blast-wave. For particles 
moving faster than the surface of the blast wave, there is a stronger tendency 
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for heavier particles to be more confined to the region of the surface owing to 
their slower thermal velocities (72). This results in heavy particles being ahead 
of lighter particles of the same asymptotic velocity, and leads to a non-zero Ar 
illustrated in Figure 3. As discussed in Section 2.3, these displacements are 
accessible through measuring odd components of the correlation function (60; 
61; 62). 

2.8 Generating Correlations Functions from Hydrodynamics and 
from Microscopic Simulations 

Any model that predicts final-state space-time and momentum information of 
emitted particles can be used to predict correlation functions. This informa- 
tion may be extracted from both microscopic simulations or from hydrodynamic 
calculations. 

Microscopic simulations model the collision by evolving particles along straight- 
line classical trajectories which are punctuated by collisions that are programmed 
to be consistent with free-space cross-sections. When the modeling is done on a 
one-to-one basis, the simulations are referred to as cascades. Boltzmann simula- 
tions are similar but employ an oversampling by a factor N s accompanied by a 
scaling down of the cross-sections by the same factor. These are then consistent 
with the Boltzmann equation and become local and relativistically covariant in 
the large N s limit (98; 99). To generate correlation functions from either class of 
simulation, there are essentially two methods which are equally justified within 
the smoothness approximation. Method I is motivated by Equation 3. This in- 
volves first creating two lists, one for each species, of the space time coordinates 
x a and Xfe of all those particles that were emitted with momenta m a P / (m a + m&) 
and mfcP/(m a + m^). From these lists, one generates <Sp(r) by sampling the 
distributions of x a — x^. This list is then convoluted with |</>(q, r)| 2 to generate 
C(P,q) for all q. In Method II, one samples pairs randomly without regard to 
their momenta. The numerator of the correlation function is then calculated by 
generating pairs with the same weight as one expects to observe experimentally 
and applying a weight given by the square of the relative wave function. The 
denominator would be calculated in a similar manner, but without the weight 
from the wave function. This method reflects the description of the correlation 
function in Equation 2. Acceptance effects or kinematic cuts can then be per- 
formed exactly as they would be performed for real particles. Method II has 
an advantage in that it is easier to accurately incorporate acceptance effects or 
tight kinematic cuts. Method I makes for a much quicker calculation because the 
procedure does not require sampling particles for irrelevant momenta (42). 

Given the equation of state and the initial energy density, hydrodynamics pro- 
vides the means for solving for the space-time development of the stress-energy 
tensor which can be used to make predictions for correlations (100; 101; 102). 
Viscous effects can also be incorporated and are non-negligible (103; 104). Gen- 
erating source functions from the output of hydrodynamic calculations is not as 
straightforward as it might seem. The Cooper-Frye prescription (105) conserves 
energy and momentum if the equation of state is one of free particles, but it 
suffers from the fact that the particles that cross backwards across the surface 
into the hydrodynamic volume enter the source function as a negative emission 
probability. If the relative velocity of the surface, as measured by an observer in 
the matter's rest frame, is not much faster than the thermal velocity, a different 
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prescription is required. Numerous prescriptions have been proposed to address 
these issues (106; 107; 108). 

Hydrodynamic models, even those that incorporate viscosity, cannot be justi- 
fied once the system expands to the point that the mean free path is similar to the 
characteristic size of the system. However, Boltzmann descriptions or cascades 
are well justified at lower densities. Several efforts have thus focused on coupling 
the two approaches (109; 110; 111; 112). because the final-state trajectories are 
established in the Boltzmann part of the prescription, one can apply either of the 
methods mentioned above. 



2.9 Phase Space Density, Entropy and Coalescence 

because phase space density depends on both the momentum p and the position 
r, a measurement of the phase space density must specify the spatial region 
over which it is determined. In practice, spatial information from two-particle 
correlation functions is instrumental to this end. For identical particle pairs, the 
"area" under the correlation function determines the average phase space density 
(113). Substituting the final phase space density for the time-integrated source 
function, 



I 



£ p (2vr)3 

and inserting into Equation 3 with Q = 2q leads to 



s(p,x)\ P0=Ep = /(p,x,t/), (26) 



/ 



d 3Q\C(Q)-l] - (2T) 3 f^|/(P/2,r, t/ )p /(P/2) 

dQ[C(Q) 1] - (2vr) [Jd3rf{p/Wf)] 2- dN/d s p > (27) 

7/ x (2tt) 3 E p dN c , 

/(p) = WTT)^ P SF ^ r =0) - 



Equation 28 applies also for the case of non-identical particles of the same phase 
space density. We note that /(p) is the phase space density averaged over coor- 
dinate space for a specific momentum using the phase space density itself as the 
weight. Unless /(p, r, f) is a constant within a fixed volume, the average phase 
space density will fall below the maximum phase space density (106). For in- 
stance, if /(p, r, t) has a Gaussian profile in coordinate space, the average phase 
space density will be 2~ 3 / 2 of the maximum phase space density for that momen- 
tum. For a Gaussian source, 

vr 3 / 2 dN 

m- W+v^L **' (28) 

where Rf nv = (E / m) -Rout-Rsidc-Riong is the product of the three radii as measured 
in the frame of the pair. The phase space density is determined by combining a 
source size measurement with the spectra. Entropy can be related to the phase 
space density in the standard way (114) 



S = (25+l)|^^[-/ln/±(l±/)ln(l±/)], (29) 

(30) 



dS/dy « J d-pr K 



dN 
d 3 p 



jj - ln(2 3 / 2 /(p)) ± /(p)/2 
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Here, Equation 30 ignores higher powers of /. 

The average phase space density is also straightforward to determine by con- 
structing ratios of spectra with species that can either bind or form a resonance. 
If species a and b can bind to form species c, thermal arguments would state 

/ c (r, t) = f a (Pm a /m c , r, t)f b (Pm b /m c , r, t)e B ' J ', (31) 

where B is the binding energy, or the excitation energy if the resonance is un- 
stable. Coalescence arguments, which give the same expression but without the 
binding energy (116; 117), are identical if the binding energy is small compared to 
the temperature, as is the case for nucleon coalescence. The average phase space 
density for a or 6 can be determined by inserting Equation 31 into Equation 28, 

t M _ - B/r m M (2g M + l) E c dN c /d 3 P 

MP) m c (2S c + l) E b , a dN b , a /d*p {6l} 

Here, the binding energy needs to be expressed in the frame of the thermal bath. 
For the case where B is small, the assumption of a thermal bath can be neglected. 
Two examples where particles with similar phase-space densities form low-energy 
resonances or bound states are pn — ► d and 4> — > K + K~ . 

By comparing the expression for / from the ratio of spectra in Equation 32, 
with Equations 28 or 28, one can determine either <Sp(r' = 0) or the Gaussian 
parameters from ratios of spectra. 

c r / = m -BjT m a m b {2S a + 1) (2S b + 1) E c dN c /d 3 P 

F[r ) (27r)3m c (25 c + 1) E a dN a /d?p a - dN/d?p b A ) 

^ (P) = (4,)^ = 0) - (M) 

Coalescence analyses can provide powerful measurements of volumes, but they 
only provide a single number, <S(r' = 0), and they cannot provide any insight 
into either the shape or the r dependence. 

3 Femtoscopic measurements 

Experimental techniques have developed considerably in response to significant 
improvements in both the theory and the quantity and quality of experimental 
data. In this section we discuss the general experimental approach for defining 
and analyzing femtoscopic correlations and their systematic dependence on global 
and kinematic quantities. 

3.1 Correlation Function Definition 

In practice, the formal definition of the correlation function in Equation 1 is 
seldom used in heavy ion physics. Instead the correlation of two particles, a and 
b for a given pair momentum P and relative momentum q, is nominally given by 
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where Ap (q) is the signal distribution, B p b (q) is the reference or background 
distribution which is ideally similar to A in all respects except for the presence 
of femtoscopic correlations, and £p(q) is a correction factor introduced to com- 
pensate for non-femtoscopic correlations present in the signal that are not fully 
accounted for in the background as well as artifacts resulting, e.g., from finite 
resolution and contamination. 

3.2 Signal Construction 

The signal Ap (q) refers to the relative momentum distribution of particles a and 
b for a given range of pair momenta, P, and a given set of event characterizations. 
Although not all analyses proceed in exactly the following fashion, the mechanics 
of constructing the signal and background are most easily understood if one 
considers as separate steps: 

1. Event quality cuts and event-class binning; 

2. Single-track (including particle identification) cuts and single-particle bin- 
ning; and 

3. two-particle pairing, two-track cuts, and pair momenta binning. 

Here the term event class refers to both physics observables, such as collision cen- 
trality and reaction plane orientation, and detector considerations, such as event 
vertex position and the condition of the detector when the event was recorded, 
usually keyed by run number. The latter considerations are relevant to the proper 
construction of the background. For event-class, particle, and pair bin, the fi- 
nal signal is usually stored as a set of 3D-histograms in the canonical relative 
momentum variables. 

Single-particle acceptances divide out with a properly constructed background, 
but 2-track acceptances can have a large effect on the correlation function. For 
this reason the analysis of 2-track cuts is the dominant consideration in the signal 
construction for most analyses. The cuts and terminology are different for Time 
Projection Chamber (TPC) experiments with near continuous hit distributions 
and Drift Chamber experiments with projective geometry, but the goals are the 
same. Split-tracks (118; 119) and ghost-tracks both refer to single tracks which 
are incorrectly reconstructed as a pair of tracks with very low relative momenta. 
Even after the event-reconstruction algorithms (which generate a list of individual 
tracks) have been optimally tuned, small traces of these false pairs remain and 
must be removed from the analysis with identical pairwise cuts. Usually only 
a tiny fraction of tracks are split, and this effect may be ignored in essentially 
all experimental analyses except femtoscopic ones. Various methods have been 
developed for identifying likely split tracks, usually based on the number (119) 
or topology (118) of space-points associated with the track. 

Pairwise effects usually also result in the loss of pairs at low relative momentum, 
because two tracks with very similar trajectories tend to be reconstructed as a 
single track. (Note that in tracking detectors, this is not a problem if one or both 
of the particles is a topologically identified neutral particle. In that case, the decay 
daughters may be well-separated even if their parents have identical momentum.) 
Such merging issues are usually resolved by pairwise cuts that remove merged 
pairs. Developing efficient and appropriate cuts can be a subtle exercise, and it 
requires good knowledge of the detector and event reconstruction software. In 
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most cases these cuts are supported by simulations, but final determinations are 
nearly always based upon data. 

What does it mean to cut out merged pairs? After all, if the tracks have merged, 
then the pair is lost anyway. The point is twofold. First, the pair efficiency 
usually does not drop from 100% to 0% sharply as a function of any variable. 
Thus, the cuts are usually tuned to exclude all (86; 118) or most (120; 121; 122) 
of the inefficient region. If regions with less than perfect efficiency remain in the 
analysis, a 2-track efficiency correction based on Monte Carlo simulations must 
be applied, typically leading to systematic uncertainties of a few percent. The 
second reason for the cut is that it is applied equally to the signal and to the 
background distribution Bp(q) (123). Thus, if some fraction / of pairs is lost 
at some relative momentum q in Ap (q), the same fraction is lost in Bp(q), and 
the ratio in Equation 35 is robust against the effect. 

3.3 Background Construction 

For reasons described above, all cut-imposed effects on the signal pair distribution 
A ab must be applied to B ab . This often means identifying which pairs would have 
been removed by merging, splitting, or other cuts, had the particles come from 
in the same event. 

The ideal background should be identical to the signal in all respects except 
for the presence of femtoscopic correlations. Therefore, the global event charac- 
teristics, single particle distributions, and acceptances should match those of the 
signal. A simple and straightforward way to construct such a background is to 
form pairs from different events within a single event class. This event-mixing 
technique (124) has gained wide acceptance in relativistic heavy-ion collisions 
where violations to energy-momentum conservation are negligible in the high 
multiplicity environment. This technique will be described in detail in what fol- 
lows. However, other methods have also been used, especially if one considers 
femtoscopy in other systems. 

For elementary-particle collisions or in low-multiplicity events, event mixing 
can violate total energy-momentum conservation, especially when exclusive fi- 
nal states or jet-axes must be preserved; thus, the correlation function would 
reflect non-femtoscopic in addition to femtoscopic correlations. In these cases, 
the most common techniques form a background from unlike-signed pairs, with 
resonance regions excluded with cuts (125) or normalized with a correlation of 
like- to unlike-signed pairs from a Monte Carlo (126). Other experiments have 
constructed a background using only Monte Carlo generated pairs (127). A 
few experiments have investigated backgrounds formed by swapping (128) or re- 
versing momentum components relative to a jet-axis (125), but these methods 
are not widely used. For detectors with symmetrical acceptance, such as the 
STAR TPC (129), momentum conservation effects may be eliminated by mixing 
pairs from the same event, with the lab momentum of one particle flipped (130). 
Backgrounds constructed from single-particle distributions as formally defined by 
Equation 1 have been used for heavy ion collisions at lower energies and shown 
to be consistent with the more commonly used event-mixing technique (131). 

In order to avoid inducing artificial structure in the correlation function, the 
particles forming pairs in the background distribution should originate from par- 
ent events with the same event characteristics. The parent events should have 
similar vertex positions to within the experimental resolution. Because detector 
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acceptances can vary with time (e.g., components may fail for some runs), parent 
events should have been measured close in time to each other; this is usually 
easiest in any case, because event mixing is done on the fly as time-ordered data 
is read sequentially. 

Parent events whose particles are mixed should also have the same single- 
particle momentum distributions. Thus, they should have similar centralities and 
orientations of the reaction plane. For example, mixing particles from events with 
very different pr slopes or directions of preferred emission (elliptic flow) would 
produce differences between A ab (q) and B ab (q) even in the absence of physical 
correlations, because almost all analyses to date have ignored these potential 
biases, it is comforting that they make little difference in practice (118). 

The list of event classes given here is by no means exhaustive. One can expect 
future analyses to incorporate the orientation of high pt particles (jet axis) or 
any other event-related observable. 

The procedure for deciding how many events to mix remains something of an 
art and involves optimizing over the range of data runs, bin width, and statistics. 
In order to minimize statistical errors, one typically forms approximately ten 
times the number of pairs in the background as in the signal. For the special case 
when all possible combinations are formed, the variance of a particular relative 

3 

momentum bin is proportional to n?, where n is the number of entries in the 
bin (22). However, as the number of pairs formed is reduced, the variance per 
bin approaches the ri2 value expected for Poisson statistics (121). It is possible 
that non- Poisson fluctuations persist in the co- variance between different bins, 
but this has not yet been investigated. 

Once the pairs have been mixed, the background must be subject to the same 
2-track cuts that have been applied to the signal. For example, the exact same 
track merging cuts or minimum separation on a detector must be applied to both 
signal and background. 

3.4 Corrections 

Corrections to the correlation function fall into three categories: finite resolution 
effects, mis-identified particle contamination, and compensation for deficiencies 
in the background. 

The first category concerns single-track momentum resolution, and reaction- 
plane resolution. We consider finite momentum resolution corrections first. Typ- 
ically, momentum resolutions are on the order of 1%. One approach is to correct 
for momentum resolution by a double ratio of the ideal correlation function gen- 
erated from a Monte Carlo simulation with perfect momentum resolution divided 
by a Monte Carlo correlation function with momentum resolution turned on. The 
femtoscopic weights are inserted into the simulations iteratively until the fitted 
radii converge (123; 119; 86; 118). A second approach is similar, but corrects 
only the Coulomb interaction term, which is most greatly affected by momentum 
resolution effects (122). In both cases, the corrections change the fitted radii by 
only <5%. 

As discussed in Section 4.2, azimuthally sensitive analyses (132; 133; 134) 
measure oscillations in correlations as a function of emission angle with respect 
to the reaction plane. Finite resolution effects in the reaction plane angle (135) 
artificially reduce the oscillation strengths. Methods have been developed to 
correct the distributions ^^(q) and B ab (q) for these effects (68; 136). 
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Another type of correction accounts for the inclusion of mis-identified and 
secondary-particle contamination. For example, electrons may be mistakenly 
identified as ir~ mesons. It is usually assumed that the mis-identified or sec- 
ondary particles are uncorrelated with other particles, so the net effect on the 
correlation function is to damp all structure uniformly in q. For purely Gaus- 
sian correlations (see Section 3.5), this effect is absorbed wholly into the A factor 
discussed in Section 2.6; homogeneity lengths — derived from the width of the 
Gaussian correlation — are unaffected by the reduction in its strength. In many 
cases, however, the homogeneity length is extracted largely from the strength of 
the correlation, and so contamination effects must be removed. In the general 
case, for which the purity p depends on the relative momentum, the correlation 
function is corrected according to C truc (q) = (C raw (q) — l)/p(q) + I (62; 118). 

It is more difficult to correct for correlated contamination. For example, if cuts 
cannot completely distinguish primary protons from those coming from A decay, 
then measured p — A correlations will contain contributions from A — A corre- 
lations. Unlike the white-noise contamination discussed above, this introduces 
structure into C pA (q) that can be accounted for only with detailed simulations. 
Such corrections will become more important at RHIC due to copious resonance 
production, and especially for baryon correlation measurements, in which the 
heavy daughter carries most of the momentum of the parent resonance. 

The last category of corrections are applied to fix deficiencies in the back- 
ground distribution. This includes corrections to account for two-particle in- 
efficiencies, which have been discussed in the previous Section. A second cor- 
rection of this type deals with the residual signal correlation that is present in 
all backgrounds derived from events that contribute to the signal. The residual 
correlation arises because femtoscopic correlations can modify the single particle 
distributions. This is especially true for small-aperture spectrometers. This effect 
can be removed with an iterative procedure (22; 137), however, for many large 
experiments the induced error is often 1% or less, and it is easier to fold this into 
the systematic errors (122). 

3.5 Fitting 

After the application of all cuts and corrections, the correlation is formed accord- 
ing to Equation 35 and then fit to determine spatial parameters. As described in 
Section 2, there are three approaches to fitting the correlation function: fitting 
to a simplified Gaussian form with strong and Coulomb interactions neglected or 
factored out, fitting to a convolution of the full kernel convoluted to a parame- 
terization of S(r), and inverting the kernel to fit a source image. The simplified 
Gaussian fits based on Equation 19 are limited to correlations of identical pions, 
kaons, and photons, but it has been the most widely used method to date be- 
cause of the computational demands of the other methods. We expect its use 
to continue for the large systematic studies in which binning in centrality, reac- 
tion plane, and kr leads to fits of more than one hundred separate correlation 
functions for a single colliding system. 

However, the functional form of the Gaussian parameterization used by exper- 
imentalists has evolved over the years. Before reviewing the most recent func- 
tional forms, it is necessary to review the treatment of the Coulomb interaction 
and the fraction of pairs coming from the source that contribute to the femto- 
scopic correlations. Both were first introduced into the literature by Zajc (22) in 
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the form of the Gamow factor given in Section 2 and an empirical parameter A 
to account for the observation that not all pairs exhibit femtoscopic correlations. 
With steady improvements in data quality and CPU speed, the Gamow factor 
has been replaced with a calculation of an squared unsymmetrized Coulomb wave 
for a finite Gaussian source. The improvements in data quality have also led to a 
self-consistent treatment of A with respect to both Coulomb and Gaussian com- 
ponents of the fit function (138; 139). For this to be accurate, we must assume 
that the source is fully chaotic, an assumption that has recently been verified 
with three-pion correlations (85; 81). The non-femtoscopic pairs consist of mis- 
identified particles and particles that emanate from too far from the source for 
the correlation to be resolved experimentally The region far from the source has 
been referred to the source halo, to differentiate it from the core. The correlation 
fit function is therefore given by Equation 36, 

C(q) = JV[AG(q)F(q) + (l-A)], (36) 

where N is the overall normalization, F is the Coulomb component, and G is the 
Gaussian form for the un-damped correlation function, Equation 19 for out-long- 
side coordinates, or Equation 21 for Yano-Koonin variables. 

Figure 4 shows projections of a 7r~ — ir~ correlation function measured by the 
STAR collaboration (118). The filled symbols are the measured correlation func- 
tion corrected for momentum resolution only and fit with Equation 36. The open 
symbols have been overcorrected by applying to all pairs the Coulomb correction 
for the fitted source dimensions. Depending on the shape of the correlation and 
degree of experimental contamination, extracted homogeneity lengths may vary 
by up to ~ 15% if the correlation function is overcorrected. 

For proton-proton correlations and non- identical particle correlations, direct 
fits are performed by convoluting the full kernel with a parameterized source. 
For these analyses, the paucity of statistics has been more of a limitation than 
the relatively modest demands in CPU power. The examples given in Section 4 
are all for one-dimensional analyses, but recent data from RHIC will soon be 
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analyzed in multi-dimensions. 

The ability to image the source by inverting the kernel is a relatively recent 
development, but one with very general applications. Because the source is pa- 
rameterized by a series of B-splines, it is a very general form which is sensitive 
to non-Gaussian shapes. To date, source imaging has been performed only with 
one-dimensional correlations, but like with the direct fits, a multi-dimensional 
kernel will soon be possible (141; 69). 

Non-Gaussian effects were reported in the first pion correlation measurement 
at RHIC (142). With higher statistics, the STAR Collaboration has studied the 
issue in greater detail (118), performing a functional expansion (the so-called 
Edgeworth expansion (143)) about a Gaussian shape. Although significant non- 
Gaussian contributions were reported, the dominant length scales were already 
extracted in the purely Gaussian fits. 

3.5.1 Minimization A simple chi-squared test is inappropriate for fitting 
correlation functions because the ratio of two Poisson distributions is not itself 
Poisson distributed, especially when taking the ratio of small numbers. For this 
reason, a log-likelihood fit function of the form given in Equation 37 is preferred. 



V 2 

X PML 



■ ( C{A + B) \ ( A + B y 



(37) 



where A, B, and Cwere introduced in Equation 35. This equation derived from 
the principle of maximum likelihood assuming that both signal and background 
are Poisson distributed (121). The full derivation of Equation 37 and comparison 
to earlier log- likelihood functions is given in (121). 



4 Measured Femtoscopic Systematics 

The first systematic study to compare femtoscopic measurements across several 
systems and experiments was performed almost 20 years ago with data from 
intermediate-energy heavy ion collisions at the Bevalac (144). The data, taken 
from experiments with different acceptances, triggers, and analysis techniques, 
were sufficient to demonstrate a crude A 1 /^ scaling of the one-dimensional radii, 
indicating that spatial scales were indeed being probed. The first femtoscopic 
measurements for relativistic heavy ion collisions were presented by the NA35 
Collaboration at the Quark Matter meeting in Nordkirchen (145; 146). More 
detailed measurements followed with the availability of sulphur and silicon beams 
at the SPS (123) and AGS (147; 148). 

Since then, increasingly sophisticated experiments at the AGS, SPS, and RHIC 
have performed femtoscopic measurements corresponding to a wide range of con- 
trol parameters. The experimental community performing the measurements has 
reached critical mass and matured substantially; a common language and knowl- 
edge base has developed concerning sometimes subtle details in performing and 
interpreting femtoscopic measurements. The result of this effort is a striking 
degree of consistency across experiments in regions of phase space where ac- 
ceptances overlap and meaningful generation of systematics across experiments. 
Large-statistics data sets routinely allow three-dimensional correlation measure- 
ments with small statistical error bars. Systematic errors, which now dominate 
the experimental errors, have been reduced to the level of ~ 5%, or ~ 0.25 fm 
for most measurements. It is no exaggeration to state that femtoscopic measure- 
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Figure 5: Pion source radius dependence on number of participants (left) and on 
charged particle multiplicity (right). Data are for for Au+Au (Pb+Pb) collisions 
at several values of y / SNN ! and also for Si+A collisions at the lowest energy. 
Average transverse momentum (&t) ~ 450 MeV/c for the PHENIX data and 
~ 390 MeV/c for the others. Data from (121; 149; 122; 118). 



ments have become a precision tool. 

Here, we cover the most important systematics of femtoscopic measurements 
from the AGS, SPS, and RHIC. We discuss only generally the physics probed 
by a given systematic, appealing to intuitive schematic models such as the blast 
wave (72). Full interpretations and comparisons to dynamic models are given in 
Section 5 

4.1 System size: iVp art and Multiplicity 

As discussed earlier, femtoscopic radii probe homogeneity regions, and not the 
entire source (hereafter, the term source will be used to refer to the entire source 
of particle emission). Nevertheless, the claim that two-particle correlations probe 
spatial scales would be given little credence if the radii did not exhibit a strong, 
positive correlation with system size. Therefore, measuring the systematic varia- 
tion of the radii vs. system composition and centrality represents the most basic 
test of both theoretical and experimental femtoscopic techniques. 

Coalescence studies (150) and two-proton measurements at the AGS (151) and 
SPS (152) unambiguously demonstrate that nucleon homogeneity lengths increase 
with decreasing impact parameter and/or increasing projectile mass, continuing 
the trend mapped at lower energies (153; 154), where directional cuts have allowed 
measurement of the shape of the homogeneity region (155; 156; 157). More de- 
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tailed information comes from pion correlations at relativistic energies, for which 
three-dimensional analyses allow partial isolation of purely geometrical effects. 
The centrality dependence of Bertsch-Pratt source radii are shown in Figure 5 
for a wide range of collision energies. The left panels show the dependence on the 
number of participating nucleons, Ap art , a generalization of the A 1 / 3 linear scal- 
ing of nuclear radii used to approximate the initial overlap geometry. All of the 

1/3 

radii exhibit a linear scaling in N p ^ rt , most with finite intercepts. Only the slope 
of the i?iong dependence shows a significant increase from the AGS to RHIC, con- 
sistent with a lifetime that increases with both centrality and ^/s^n. The trend 
of increasing i?i ong with increasing ^/snn is reversed for ^/snn < 5 GeV (86) . 

The right panels of Figure 5 show the same radii as a function of (dN^/drj) 1 ^ 3 . 
The primary motivation for exploring the (ciiV c h/(i?7) 1 / 3 dependence is its relation 
to the final state geometry through the density at freeze-out. However, the two 
scaling quantities are highly correlated. In fact, the values of dN c h/dr] shown 
on the right side of Figure 5 were derived from iVp art using the N® t parame- 
terizations given in (158), and conversely, the N part values are often calculated 
from multiplicity distributions using a Glauber model. Given this caveat, the 
inside and -Riong values exhibit a linear dependence on (dN^/drf) 1 ^ , again with 
finite intercepts. The strong uniformity from ^/snn °f 5 to 200 GeV leads one to 
believe that the approximate -/Vp art scaling (initial overlap geometry) is a result 
of the scaling with multiplicity (final freeze-out geometry) and not the other way 
around. 

The parameter R ou t , which mixes spatial and temporal information (see Section 
2.5), increases with multiplicity at each given collision energy, but does not follow 
a universal curve. However, the strikingly ^/.SNN-independent multiplicity scaling 
of the geometric radii -R s idc an d i?i on gStrongly suggests that the observed increase 
of these radii with collision energy for ^/snn > 5 GeV (see Section 5.1) is due 
simply to the rise of multiplicity with collision energy. This trend, as well as its 
violation at ^/snn < 5 GeV, has been interpreted in terms of changing chemical 
composition of the source as the system evolves with energy from baryon to meson 
dominance (210). 

We note that the systematics in system size represent an initial sanity check 
for the femtoscopic technique. The obvious direct connection of the radii to the 
source geometry estimated in two ways refutes suggestions (159) that smaller, 
non-geometric length scales dominate experimentally extracted transverse radii. 

4.2 Source Shape: Pair emission angle relative to b 

The variation of femtoscopic radii with the pair emission angle relative to b (</> pa i r ) 
can be used to probe the three-dimensional shape of the source (75; 160; 66; 161; 
162; 67; 68; 72). The anisotropic shape transverse to the beam direction — the 
coordinate-space analog to the elliptic flow characterizing momentum-space — 
gives rise to cos(n</> pa i r ) (n even) oscillations in the squared transverse source 
radii R 2 out , R 2 sidc , inside (66; 68). 

Just as one expects the source (and homogeneity regions) to be larger for 
decreasing \b\, one also expects it to be rounder, reflected by small oscillations 
of the radii. Figure 6 for mid-rapidity pions from Au+Au collisions at RHIC 
confirms this expectation. As |6| increases, the oscillations indicate a transverse 
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Figure 6: Squared source radii mea- 
sured at mid-rapidity in Au+Au col- 
lisions at RHIC, as a function of pair 
emission angle relative to the reaction 
plane. Data for three centralities are 
shown. Figure from (134). 



Figure 7: Squared source radii mea- 
sured near mid-rapidity in mid-central 
Au+Au collisions at lowest-energy AGS 
settings, as a function of pair emission 
angle relative to the reaction plane. Fig- 
ure from (132). 



source increasingly elongated out of the reaction plane (134). 1 

The strong in-plane expansion (163) does not fully convert the initial out-of- 
plane (overlap) geometry into an in-plane-extended source at freeze-out. This 
suggests a rather short evolution time; in essence, the system did not have time 
to reverse its deformation. However, this is only a hint, and a full dynamical 
transport calculation is required to extract physical timescales (90). 

Whereas at the highest RHIC energy, the freeze-out anisotropy is ~ ^ of the 
initial (134), at low AGS energies, the final anisotropy is consistent with that of 
the initial overlap region (132), or perhaps slightly lower, because elliptic flow 
vanishes — changes sign — at these energies (164), these trends make intuitive sense 
and suggest an underlying connection to the evolution dynamics. It would be 
desirable to map the source anisotropy at intermediate (AGS and SPS) energies, 
for which there may be interesting changes in the space-time systematics. At 
these energies, there have been intriguing hints of asymmetries in the homogeneity 
regions for pions (165; 166; 167; 168) and protons (169), and in the proton-pion 
separation (166; 73), although they have not been finalized. 

If the impact parameter direction b — not simply the 2 nd -order event-plane angle 
(unambiguous only over a range [0, tt]) — is known, then more detailed information 
may be obtained. In the left panel of Figure 2, the source is tilted with respect to 
the beam axis, toward b. Just as anisotropic azimuthal geometry in the transverse 
plane is related to the structure of elliptic flow (90; 100), a tilted geometry can 
reveal important information on the underlying nature of directed flow (170; 

lr rhe out-of-plane nature of the elongation may be read directly from Figure 6. Ignor- 
ing collective flow or opacity effects (e.g. 67) an out-of-plane-extended Source would produce 
^sideWpair = 0°) > R^idc (^pair = 90°), as seen in Figure 6. Collective flow effects complicate 
this picture (66; 72), but the sign of the oscillations are determined by geometric, not dynamic, 
effects for realistic sources at RHIC (90; 72). 
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Figure 8: The Yano- 
Koonin source rapidity 
is shown as a func- 
tion of the pion pair 
rapidity for central 
Au(Pb)+Au(Pb) col- 
lisions over a broad 
range of energies (open 
symbols are for ir~ir~, 
closed symbols for 
7r + 7r + ). Both quantities 
are in the center-of- 
mass frame of the 
colliding system. Data 
from (176; 177; 178; 
179). 



2 - 



1 - 





nJI 


<k T > 




■ 


5 


0.10 GeV 


E877 





17 


0.15 GeV 


NA49 




17 


0.45 GeV 


WA97 


• 


62 


0.45 GeV 


PHOBOS 




200 


0.45 GeV 


PHOBOS 



A 



171; 67; 132). The structure in R 2 ol and R 2 sl shown in Figure 7 is generated by 
this tilt. The spatial tilt has been measured only at low AGS energies (132); a 
measurement at RHIC might reveal exotic geometric configurations generated by 
quark gluon plasma (QGP) formation (171; 172) and would impact the important 
issue of boost-invariance at mid-rapidity. 

4.3 Boost invariance : Y T7T 

In high-energy hadronic collisions, the initial parton distribution is expected to 
be approximately flat in rapidity. This momentum rapidity distribution may 
correspond to producing matter that initially exhibits a boost-invariant Hubble- 
type scaling correlation between longitudinal flow velocity and space-time points, 
vljiow = z/t (173; 174). Relativistic hydrodynamics preserves boost-invariance 
of the initial conditions throughout its dynamical evolution (174). The combi- 
nation of the above arguments underlies expectations that in ultra-relativistic 
heavy ion collisions, particle production emerges from boost-invariant longitu- 
dinal flow, and that dN/dy exhibits an approximately boost-invariant plateau 
around mid-rapidity. However, an extended plateau has never been observed 
from AGS through RHIC energies (3). 

because correlation measurements access spatio-temporal information, the ques- 
tion arises (78; 175) whether they allow us to test the relation vlji ow = z/t 
between the space-time rapidity and the momentum rapidity of the source. For 
boost-invariant sources, one can show that the pair momentum rapidity is equal 
to the Yano-Koonin source velocity, which is directly obtained from the Gaussian 
radius parameters, Y n7T Yljiow = tanh -1 vljiow- However, even if the source 
density distribution shows significant deviations from boost-invariance, this re- 
lation still holds approximately as long as the velocity profile is boost invariant, 
and kx is sufficiently large (78). 

Figure 8 reveals a roughly universal dependence of Yyk on Y W7r for pions 
from central collisions, depending weakly, if at all, on ^/snn (see Section 2.5.3). 
This trend is particularly striking given the very different center-of-mass projec- 
tile rapidities (~ 1.55 and 5.5 for ^/snn = 5 and 200 GeV, respectively) and 
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corresponding widths of the pion distributions dN/dy. 

By way of caution, we note that the results from RHIC are limited to the 
region Y W7T < 1.2, and that the deviations from boost-invariance are mostly in 
the lower energy data. Extending the RHIC results to more forward rapidities 
would provide a important test for both the velocity scaling at RHIC and the 
energy-independence that is exhibited in Figure 8. 

For central collisions, the roughly universal behavior approximately obeys the 
boost-invariant consistency relationship discussed above. Moreover, Yljiow shows 
a significant kx dependence and falls below the linear relation in particular for 
small kx (178). Qualitatively, this is consistent with blast-wave models in which 
a boost-invariant longitudinal flow is superimposed on a source density distri- 
bution of finite longitudinal width. However, a full dynamical understanding of 
the dependence is missing so far. The flat YVk dependence on Y W7T measured at 
the SPS (178) for the most peripheral collisions is counter-intuitive, and requires 
further study. 

The question of whether the source has boost-invariant space-time structure 
is an important one. There are many reports of very short evolution timescales 
(lifetimes) based on fits to the data with Equation 24, which is based upon an 
assumption of boost-invariance (87). Relaxation of that assumption might lead 
to considerably larger estimates (180). 

4.4 Collective dynamics: k T and particle mass 

As discussed in Section 2.7.2, the dynamic substructure of the source is encoded 
in space-momentum (xp) correlations. Longitudinal xp correlations, encoded in 
-Riong(^r)> are generally acknowledged (87; 72) to reflect longitudinal flow, be- 
cause all transverse correlations are generated in the collision itself, considerably 
more attention has generally been paid to the transverse substructure than to 
the longitudinal flow discussed in Section 4.3. 

The most common explanation for transverse xp correlations is collective trans- 
verse flow (96). These correlations have mostly been studied through pion corre- 
lations, but transverse flow implies also a systematic trend as the particle mass 
is varied. 

4.4.1 kx dependence OF pion radii Collective flow generates a charac- 
teristic fall-off of the pion source radii with kx, which is ubiquitously observed in 
data. Final results for the fey-dependence of Gaussian radii from central Au+Au 
(Pb+Pb) collisions exist at the AGS (86; 121), SPS (177; 137; 181; 178; 149), 
and RHIC (142; 182; 134; 118; 179; 122). As is clear from Figure 9, aside from a 
small variation in overall scale (discussed later) , the kr dependence is startlingly 
similar for all energies. 

Figure 10 quantifies the evolution of the my-dependence of the pion source 
radii with ^/snn> using fits to Ri{m,T) ~ m^"'. As discussed in Section 2.7.1, 
on = 0.5 would represent expectations for instantaneous thermal emission for a 
three-dimensionally expanding fireball in the limit of large tut- 

The similarity persists as N pait is varied. In Au+Au collisions at RHIC, 

1/3 

Ri{ m r) {i = out, side, long) simply scale as N' t , with perhaps some flattening 
for -/Vpart < 100 (118; 122). Very similar kT dependence for different iV part is also 
observed in Pb+Pb collisions at SPS (149) and for Si+Au and Au+Au collisions 
at the AGS (121). 
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Figure 9: World data set of published dependence of pion Bertsch-Pratt 
radii near mid-rapidity from Au+Au (Pb+Pb) collisions. Centrality selection is 
roughly top 10% of cross-section, but varies somewhat with experiment. Data 
from (86; 121; 149; 137; 178; 181; 182; 122; 179; 142; 118) Lines represent 
parameterized fits; see text for details. 



In a flow-dominated freeze-out scenario, the fall-off of transverse radii with tut 
increases as flow increases and/or temperature decreases (e.g. 72). Blast-wave fits 
to spectra (183) indicate that freeze-out flow and temperature vary significantly 
with y / SNN for yj snn % 10 GeV. The overall approximate -^/sNN-independence 
of the aj parameters may reflect the fact that significantly changing the slope of 
Ri{rriT) requires very large changes in flow and temperature; on the other hand, 
it could be that the compensating effects of smaller (larger) homogeneity lengths 
generated by larger flow (temperature) cancel almost exactly in nature. Although 
Ri{rriT) almost certainly reflects strong collective flow, determining the strength 
of that flow requires other information, such as particle spectra (94; 72). 

because the radii fall off roughly as 1/ \frrvr (see Figure 10) and such a de- 
pendence has been discussed frequently in the literature (e.g. 184; 185), it is 
interesting to look at the overall scale parameter from a single-parameter fit to 
R\j y/mT- The y / SNN-dependence of -Rg ide and -R( ong are shown in Figure 11. The 
scale of the longitudinal homogeneity length grows significantly with y / SNN, con- 
sistent with an increase of the system evolution time. However, R' varies only 
very weakly with ^/snn- 

4.4.2 Systematics with particle mass Systematic studies for different 
mass particles provide additional controls probing the space-time evolution of the 
source. In particular for kaons, the interpretation may be simplified owing to re- 
duced effects of resonance feed-down (186) and a reduced scattering cross-section 
for K + in nuclear matter, raising the possibility that kaon correlations could peer 
farther back to earlier stages of the collision (187). Indeed, the first kaon mea- 
surements (188; 189; 190; 176) reported smaller source radii for kaons. However, 
the observation that radii for K + and K~ were very similar (189) was an early 
experimental indication that different cross-sections were not the driving physics 
behind these smaller radii. This was supported by model calculations (191), 
which suggested that K + and K~ in fact scattered roughly equally in the dense 
medium created in heavy ion collisions. In this case, the smaller radii for kaons 
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Figure 10: The excitation function of 
nix power-law fall-off of pion source 
radii (i^(mx) ~ m^ a '). Shaded regions 
show cti = 0.5 ± 0.05. 



Figure 11: The excitation function of 
scale parameters R[ from fits of the 
data in Figure 9 to iJj(mr) = R'Jy/rfvr. 
Note: this is a different functional form 
than that used in Figure 10. 



results from their increased mass in a flow field, not different cross-sections. 

If indeed flow is generated in matter sufficiently dense that individual cross- 
sections are unimportant, then all particles participate equally in collective trans- 
verse flow. In this case, their source radii should approximately follow a common 
m,T scaling (30; 159; 192; 8). Within uncertainties, first results on kaon inter- 
ferometry by NA44 at the SPS in S+Pb (184) and Pb+PB (193) collisions were 
consistent with a 1 / ^Jm? scaling expected for isotropic Hubble flow (30; 159; 192). 
In some more recent analyses (194), the common mj- systematic for the trans- 
verse radii is less steep than for R\ onz {mT) (see Figure 10), as might be expected 
for more boost-invariant (non-isotropic) flowing systems (8; 72). 

Figure 12 collects the dependence of homogeneity lengths for several en- 
ergies. The left panels show results for Si+Au collisions at ^/snn=5.4 GeV, 
measured by E802 for pions (121) and kaons (188; 190). Femtoscopic radii for 
pions (137; 194; 149; 181) kaons (193), protons (195), and photons (196) mea- 
sured in Pb+Pb collisions at the SPS are shown in the center panels. The right 
panel shows the one-dimensional radius parameter i?j nv measured at RHIC for 
pions, charged kaons, and protons (197), neutral kaons (198), and with A — p 
correlations (199). 2 To compare across energies, -Ri nv results are included for the 
AGS and also for the SPS, where the i?i nv values were calculated from the 3D fit 
results by accounting for the boost along the outwards direction from the LCMS 

(in which P z = 0) to the PCOM (pair center-of-mass, in which \P\= 0) frame, 
-ftfnv = ^long + -^side + T^out) where 7 is given by m T /m of the pair. Note 
that for massless particles, such as photons, 7 is given by kx/Qmv For a given 
orientation of the photon pair momentum i?; nv can be related to 3D radii in the 

2 Protons and A baryons are not identical particles, of course. However, their masses are 
sufficiently close to try including Rf^ on a mr-scaling plot. Adherence to the scaling is 
consistent with a flow-dominated scenario in which homogeneity lengths depend only on mass. 
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Figure 12: Transverse mass dependence of homogeneity lengths from correlations 
between particles of (nearly) identical mass. 

LCMS frame through the following expression for the exponent of the Gaussian 
source term, 

QiLv^iLv = Qinv { R lon g cos2 $ + ^side sin2 Osin 2 (f> + fcf./ Q? nv ^out sin2 0cOS 2 (f) 

= k^Rl ut sin 2 6cos 2 4> + Qf nv (Rf ong cos 2 6 + R 2 ide sin 2 9sin 2 <j)) . (38) 

Thus .Rout drops completely out of the source term and enters into the A coeffi- 
cient (196), but as a function of hp. For this reason the ID Q mY fits from WA98 
are plotted alongside values of i? s ide and R\ oa g instead. 

This consistency between different particle types may carry an important mes- 
sage. It calls into question theoretical scenarios which appear to explain Ri{rriT) 
for particular particle types only (33). Further, the consistency with emission 
from a common flow-dominated source may also support freeze-out scenarios in 
which the last scattering in the dense phase determines the homogeneity region, 
instead of milder rescatterings in the more dilute stage, which are dominated by 
particle species-dependent cross-sections (200; 201). 

4.5 New possibilities: Non-identical particle correlations 

Almost all femtoscopic measurements to date have been done through correla- 
tions of identical particles, usually charged pions. With the availability of high- 
statistics data sets and new theoretical ideas (202), however, experiments are 
beginning to make full multi-dimensional analyses of non-identical particle cor- 
relations. These correlations are being used to test and refine the treatment of 
Coulomb effects in identical-particle correlation analyses, to explore violations of 
flow-induced universal tut scaling, and to extract qualitatively new information 
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on the space-time substructure of the source. 

4.5.1 Opposite-sign pion correlations The most common non-identical 
particle femtoscopy measurements have involved n + — n~ correlations. Here, the 
natural assumption is that the homogeneity regions of these particles coincide. 
The primary interest is in whether the homogeneity volume is extracted from 
^ — 7r =l= correlations — driven by quantum statistics and Coulomb final state in- 
teractions (with minimal strong FSI effects)- is consistent with that extracted 
from 7r + — 7r~ correlations — determined only by Coulomb and strong FSI (118). 
This issue has obvious implications for the FSI weight F(q) discussed in Sec- 
tion 3.5. 

At the SPS, 7r + — 7T~ correlations from reactions with sulphur (203) and 
lead (177) beams were consistent with emission from a homogeneity region with 
spatial scale roughly the same as the radii from identical-particle correlations. 
Preliminary analysis of opposite-sign pion correlations from Au+Au collisions at 
the AGS (176; 73) similarly find consistency with scales extracted from like-sign 
pion correlations. Recently, STAR reported (118) strong consistency at RHIC 
energies as well. Using large-statistics data sets, they further show that contribu- 
tions to the 7r + — 7T _ correlation function from strong-force interactions, though 
small, are nevertheless important in explaining the data in detail. 

4.5.2 Other non-identical particle correlations in \q\ For statisti- 
cal reasons, most non-identical particle correlations are measured in one-dimensional 
|gj-space, and thus probe only the average spatial separation between the two 
particles at freeze-out in the pair center-of-momentum frame (see Section 2.3). 
Because they are sensitive to the size of each particles emission source and this 
separation between them, it is often necessary to rely on the identical particle 
correlations for a proper interpretation of results. 

Despite these ambiguities, one-dimensional non-identical particle correlations 
may be used to test existing systematics and expectations of the freeze-out sce- 
nario. Conversely, if the freeze-out geometry may be taken as given, non-identical 
particle correlations may place constraints on hadronic scattering parameters, e.g. 
to measure the squared relative wave-function of a pion and a E baryon. First 
studies along these lines are underway at RHIC (204). 

Wang and Pratt (58) suggested the measurement of p — A correlations, which 
may be more sensitive to large structure than p— p correlations, and which have a 
higher two-track reconstruction efficiency for some experiments. Published results 
from Au+Au collisions at the AGS (55) suggest a p — A separation distribution 
with a width similar to the proton homogeneity length and roughly consistent 
with tut scaling expectations (cf Section 4.4.2). Preliminary results at SPS (205) 
and RHIC (204; 206) give similar conclusions, because the p — A potential is not 
unambiguously known theoretically (58), it is unclear whether possible statisti- 
cally marginal violations of scaling (55) or differences between baryon and 
anti-baryon emission regions (204) are meaningful. To first order, existing p — A 
correlations confirm existing systematics. 

With the high-statistics and high-quality data sets at the largest collision ener- 
gies, truly exotic correlation studies are possible. Preliminary results (206) from 
RHIC on 7r — S correlations look particularly promising. Here, the well-measured 
pion emission distribution may be used to study the strange baryon freeze-out 
configuration. It may also provide information on the tt — E final state interaction 
and scattering cross-section which in turn can be used to constrain our under- 
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standing of the sources of collective flow. Blast-wave calculations reproduce the 
preliminary tt — S correlations, suggesting that the S flow is determined by its 
mass, not its quark content. 

4.5.3 New information: non-identical particle correlations with 
directional CUTS ON q Non-identical particle correlation analyses as a func- 
tion of — or with cuts on — the relative direction of q and P reveal qualitatively 
new information (60; 61; 202, and Section 2). In particular, in addition to the 
root mean square (RMS) width of the separation distribution, the direction and 
size of the average separation between the particles is probed, although offsets 
in time and space cannot be disentangled; this is shown as Ar in Figure 3. The 
correlation functions selected for q || P and q P differ if Ar 7^ 0. Further- 
more, collective flow will induce position-momentum correlations detectable with 
directionally selected non-identical particle correlations (61; 72). 

These correlations are statistically challenging, and few results are available. 
At RHIC, STAR has reported (62) asymmetries in — ir^ correlations mea- 
sured in central Au+Au collisions. Blast- wave calculations with transverse flow 
roughly adjusted to reproduce other observations at RHIC (72) describe the data 
semi-quantitatively. Preliminary studies of p— tt correlations at the SPS (205) and 
RHIC (206) exhibit similar mass-ordered spatial asymmetries in the transverse 
plane. A preliminary study at the AGS (73) reported very large (Ari ong ~ 10 fm) 
average p — tt separations in the longitudinal (beam) and impact parameter direc- 
tion for forward-moving particles, suggesting strong longitudinal flow; however, 
this result was never confirmed. 

We expect full three-dimensional analyses of a wide range of non-identical par- 
ticle combinations to be available in the near future from RHIC experiments. 
Sophisticated analyses may probe non-trivial geometric substructure when se- 
lecting on reaction plane, including the sideward shift predicted by blast-wave 
calculations (72) when anisotropic flow structure is present. 

5 Interpretations of Experimental Results 

In this section, we ask what we can learn from the spectrum of results just 
presented. Beginning with the broadest, least detailed observations, we move to 
two fundamental quantities that may be directly extracted from the data, and 
finally on to comparisons with specific models of heavy ion collisions. 

5.1 General Conclusions from Systematic Trends 

One of the first messages to take away from the discussion of Section 4 is that 
the results are stable across detector and method; experimental systematics and 
uncertainties are under control. Whatever difficulties we may have in interpret- 
ing measurements, we may be confident that they do not have their origin in 
experimental artifact. 

The size and shape inferred from two-particle correlations tracked with colli- 
sion geometry as anticipated. Kinematic and mass dependences of femtoscopic 
measurements showed the expected clear signatures of strong collective flow in 
the beam direction and perpendicular to it. 

At a generic level (ignoring quantitative predictions), we are first given pause at 
the jejune nature of the ^Js n Ar-dependence of femtoscopic parameters. The most 
common example discussed is the pion source radii excitation function shown in 
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Figure 13. In Section 4 we explored the trends in considerably greater detail, but 
the figure conveys the right message: Scanning ^/sJfN through a range of two 
orders of magnitude changes final-state geometry little. 

Based on rather generic arguments of soft points in the equation of state or 
entropy generation during a phase transition, there had been hopes for non-trivial 
structure in the excitation function, as the energy threshold for quark gluon 
plasma (QGP) creation was crossed (63; 207; 208; 209). General expectations 
were for long (relative to the explosion of a purely hadronic system) system 
evolution time scales if QGP was formed. Blast-wave analyses of data (e.g. 72) 
appear to rule out systems with lifetimes in the neighborhood of 20 fm/c or 
higher. However, owing to dynamic effects, careful comparison with a dynamical 
model is required to extract detailed evolution information; this is discussed in 
Section 5.3. 

As discussed in Section 4.1, the weak increase with i? out of the femtoscopic 
radii for central collisions shown in Figure 13 is mostly due simply to increasing 
multiplicity. The CERES collaboration (210). has suggested that this increase 
may be understood in terms of an energy-independent mean free path of 1 fm at 
freeze-out. Interestingly, this explanation also appears to describe the decrease 
of radii with R out < 5GeV, when the different cross-sections and abundances of 
pions and protons are accounted for. This observation emphasizes the importance 
of chemistry in the determination of the freeze-out geometry. However, it neglects 
the dynamical structure of the source (e.g., flow) and the importance of the 
six-dimensional phase space density discussed in Section 2.9. We discuss this 
fundamental quantity next. 

5.2 Phase Space Density and Entropy 

As shown in Section 2.9, average phase space densities, /(p), can be calculated 
by combining source-size measurements with spectra. Pionic phase space den- 
sities have been estimated for 130^4 GeV collisions at RHIC (114), at the top 
SPS energies, and for several AGS energies (213). Figure 14 shows results for all 
three regions. In each case, the phase space densities were calculated via Equa- 
tion 28. For the SPS case, results were generated by applying Equation 28 to 
published spectra (214) and source-size measurements (215). We note that our 
calculations for the SPS are higher than previously published values at low pt 
(216). This discrepancy is likely due to the fact that in Reference (216), analytic 
parameterizations were used that differ significantly from published spectra at 
low pt ■ 

When applying Equation 28, an issue arises as to whether one should subtract 
the contribution from resonances to the spectra. Indeed, if the pions are created 
by decays so far outside the source volume that they do not contribute to the 
correlation function, they should not be considered as pions when calculating the 
phase space density. There are two strategies to correct for such pions. First, 
one could use spectra where such particles are subtracted and apply Equation 28 
literally. As a second option, one could use the spectra without subtractions, but 
then multiply the expression for / by \f\. because most published spectra have 
been purged of the products of weak decays, the first method is usually applied. 
However, published spectra still include the contribution from 77's, which decay 
thousands of fm away from the source. The 77 contribution was accounted for in 
the SPS calculation for Figure 14 by reducing the spectra by 5%. 
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Figure 13: Panels (a-d): 
Excitation function of 
7r~ source parameters 
at mid-rapidity and low 
k T (~ 0.17 GeV/c) in 
central Au+Au(Pb+Pb) 
collisions. PHENIX data 
are for hp ~ 0.26 GeV/c 
and so fall somewhat 
below the trend. Panel 
(e): Transverse freeze-out 
anisotropy parame- 
ter from non-central 
(|6| ~ 8 fm) Au+Au 
collisions, estimated from 
the azimuthal depen- 
dence of source radii. 
Data are from (182; 122; 
142; 134; 118; 121; 149; 
181; 179; 137; 86; 132). 
Also shown are calcula- 
tions (86; 132; 137; 211) 
at several energies with 
the RQMD model (212); 
hashed region at other 
values of ^snn inter- 
polates between these 
calculations. 
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Figure 14: Aver- 
age pionic phase 
space densities for 
central Au+Au and 
Pb+Pb collisions 
from the AGS to 
RHIC rise with 
beam energy but 
seem to plateau at 
SPS energies. Val- 
ues were calculated 
with Equation 28. 



0.4 



0.2 - 



0.0 




0.0 



0.2 0.4 
P T (GeV/c) 



The phase-space densities in Figure 14 show a steady rise with beam energy 
that seems to plateau at SPS energies, because the displayed phase space densities 
have been averaged over coordinate space, the peak values are higher (91; 106), 
by a factor of 2\/2, if the spatial profiles are Gaussian. For a breakup temperature 
of 110 MeV, this requires a rather high chemical potential, near or above 80 MeV 
(217; 218; 219). 
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Ratios of particle yields are consistent with chemical freeze-out at temperatures 
near 170 MeV (220; 221; 222). This suggests an interpretation where at higher 
densities the system is so strongly interacting that yields equilibrate until the 
system reaches this temperature and are then frozen during the subsequent freeze- 
out. At AGS energies, the increase in pion production brought on by increased 
beam energy results in more pions being pushed into a given amount of phase 
space as realized in Figure 14 by the rising phase space density. However, one 
would expect that as the excitation energy surpasses the threshold for a 170 
MeV temperature, the phase space density would reach a limiting value, because 
the average phase space density is preserved in an isentropic expansion with fixed 
particle number (not exactly true if several masses are mixed together) , one would 
expect phase space densities to saturate once energy densities reached this value. 
Indeed, the behavior in Figure 14 is consistent with this scenario. 

As shown in Section 2.9, entropy can be calculated from average phase space 
densities using Equation 30. Spectra and source-size measurements for baryons 
and mesons were used to estimate phase space densities and entropy for 130A 
GeV Au+Au collisions at RHIC. The total entropy in the central unit of rapidity 
was estimated at dS/dy = 4450 ± 10% (114). In a hydrodynamic expansion, 
entropy is conserved though viscosity, and shock waves might result in a roughly 
10% increase during the evolution. Thus, this measurement provides an upper 
bound for the entropy at r ~ 1 fm/c when thermalization first occurs. At r = 1 
fm/c, the volume for particles in this rapidity slice is determined by the geometric 
cross-sectional areas of the overlapping gold nuclei multiplied by cr, thus an upper 
bound for the entropy density can be determined s < {dS / dy) / {tit R 2 ) . Knowing 
the energy density at r = 1 fm/c would then provide a point in the equation 
of state, s vs. e. The value of dS/dy estimated in (114) is consistent with 
lattice calculations if e(r = 1 fm/c) ~ 7 GeV/fm 3 . Estimates of the original 
energy density from the final-state measurement of dEt/dn are in the range of 
4.5 GeV/fm 3 (223), but because these estimates neglect losses from longitudinal 
work or the energy from longitudinal thermal motion, the 7 GeV/fm 3 value seems 
reasonable. 

5.3 Dynamic Models and their Comparison with Data 

It is increasingly recognized that the comparison of dynamic models of heavy ion 
collisions to data is only insightful if it involves a sufficiently large variety of ex- 
perimental data. A comparison of dynamic models to femtoscopic measurements 
alone (or to any other class of measurements alone) is of limited value, simply 
because for realistic models, the number of possible model-dependent parameter 
choices then tends to exceed the number of experimental constraints. In fact, all 
the model results that we review in the current subsection remain unsatisfactory 
with this respect: They either deviate significantly from femtoscopic data, or 
they reproduce these data at the price of missing other important experimental 
information. In particular, there is so far no dynamically consistent model that 
reproduces quantitatively both the systematic trends discussed in Section 4 and 
the corresponding single inclusive spectra. In this situation, the scope of this sub- 
section is somewhat limited. We want to explain why a dynamical understanding 
of femtoscopic measurements is important. We shall also discuss the key physics 
input that enters current attempts of dynamic modeling and the uncertainties 
resulting from it. However, we shall try to bypass as far as possible model- 



38 



Lisa et al. 



dependent details and rather focus on the question of which qualitative changes 
in the underlying dynamics result in characteristic changes of femtoscopic data. 

Correlation measurements provide a snapshot of the geometrical distribution of 
particles at the time they decouple from the reaction. This geometrical distribu- 
tion provides a unique test of the dynamical evolution of the produced matter at 
the late stage, because the spatial extension, dynamical evolution, and lifetime of 
the produced system determines phase space density and thus particle reaction 
rates, any dynamic model for the latter has to be consistent with femtoscopic 
information. So far, correlation analyses have focused mainly on Boltzmann (or 
cascade) models, on hydrodynamic models, or on combinations of both (hybrid 
models). These model classes correspond to rather different equations of state. 

The equations of state represented by cascade models tend to be stiff unless 
they incorporate large number of resonant scatterings. If particles collide via 
2 — > 1 — ► 2 processes where the intermediate state has a finite lifetime, the equa- 
tion of state can be softened (224; 225). A prominent example of a cascade model 
is RQMD (191), Relativistic Quantum Molecular Dynamics; it is the only one 
which has been compared to data at AGS, SPS, and RHIC. Results for RQMD 
are shown in Figure 13. However, for RHIC energies, the dynamic consistency of 
RQMD is questionable because for most of its evolution, the model uses hadronic 
degrees of freedom, although RQMD simulations for RHIC yield an energy den- 
sity which stays above that of normal nuclear matter for a significant duration 
(~ 5 fm/c). As can be seen in Figure 15, RQMD, which models the expan- 
sion as a hadronic cascade, overpredicts i? ou t and i?i on g at RHIC despite the fact 
that it underpredicts multiplicities. Another hadron cascade model, the Hadronic 
Rescattering Model (226; 227), based solely on hadronic rescattering with sudden 
collisions, gives smaller sources, which come closer to the data. Part of the dif- 
ference between the two hadronic cascades may derive from RQMD's treatment 
of scatterings as resonant interactions with finite lifetimes, which differs from the 
instantaneous collisions employed in the Hadronic Scattering Model. Figure 15 
also displays results for Molnar's Parton Cascase (MPC) (228), which aims to 
provide a transparent partonic toy model by modeling the collision of light par- 
tons which undergo a one-to-one hadronization to pions. MPC, which has only 
instantaneous two-to-two scatterings, should have the stiffest effective equation 
of state and underestimates the source radii. 

On of the most complicated transport codes available is A Multi-Phase Trans- 
port (AMPT) model. AMPT (229; 230) aims at a realistic description of all 
aspects of the reaction dynamics and includes a partonic cascade coupled to a 
hadronic cascade employing a full list of resonant interactions. AMPT provides a 
good fit to experimental radii, though the kr fall-off is stronger than that found 

— 1/2 

in the data, which fall off ~ m T . The more rapid decrease of radius parameters 
with respect to kx may be due to the continuous surface-like emission charac- 
teristic of microscopic models, which when combined with cooling gives a higher 
relative weight for high-energy particles to have been emitted earlier in time, 
before the reaction volume has reached its full spatial extent. Any modification 
that would reduce this type of emission should improve agreement with data. 
Cautiously, one would conclude that the results in Figure 15 favor models with 
a stiff, but not too stiff, equation of state and no latent heat. 

This conclusion is further supported if one compares the results of cascade 
models with those of hydrodynamic simulations or hybrid models. Figure 16 
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Figure 15: Pion radius parameters from four Boltzmann/cascade models are 
compared to experimental RHIC data. For Molnar's Parton Cascase calculations, 
the reported radii are simply space-time variances. For the other models, radii 
are obtained from Gaussian fits to correlation functions generated according to 
the methods described in Section 2.8. Measured data for positive and negative 
pions are indicated by closed and open symbols, respectively. 

compares results from RHIC to a three-dimensional hydrodynamic model of Hi- 
rano et al. (231) that investigates the effects of resonance decays on chemical 
composition, the two-dimensional model of Heinz and Kolb (232), and a two- 
dimensional chiral model by Zseische et al. (233) that performs calculations for 
both first order and cross-over transitions. In all cases, the more favorable calcu- 
lations are compared with the data: partial chemical freeze-out for Hirano, and 
cross-over transition for Zschiesche, but large discrepancies between the mod- 
els and the data still remain. The results from Zseische were taken from the 
cross-over transition for the lowest critical temperature, 80 MeV. Substituting 
a first order phase transition or increasing the critical temperature to 130 MeV 
increased the value of -Rout relative to -Rsidej thereby increasing the differences 
with the data. Hirano also reported a calculation under complete chemical equi- 
librium at freeze-out; this improves the agreement with V2, but at the expense of 
a much larger value of -R ou t/-R S ide- 

In general, these models invoke equations of state that are typically softer than 
those used in cascades and Boltzmann calculations, and they often have latent 
heats to accompany the transition from the partonic phase. As a consequence, 
lifetime and emission duration of the produced matter are significantly larger 
than what one finds in cascade models, and such models often significantly over- 
predict i?ion g and -Rout/-Rside- The fact that R sl de comes out smaller than the 
data in Figure 16 is mainly due to an attempt to compensate within the available 
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Figure 16: Hydrodynamic (Zschiesche, Hirano, and Kolb) and hybrid hydrody- 
namic/cascade (Soft 7 ) models calculations in comparison to RHIC data. Data are 
for 27r - (open symbols) and 2tt + (closed symbols) source radii. 



model-parameter range the very large time-scales as much as possible. Purely 
hydrodynamic models can vary overall source volumes by adjusting the break-up 
criteria, but doing so can make it difficult to fit the three source dimensions, 
and their txit dependence, as well as other observables. Several prescriptions 
have been applied to improve the modeling of the breakup in the late stage 
to depend on microscopic considerations determined by free-space cross-sections 
without having an extra adjustable parameter, e.g. break-up density, (107; 108; 
106). An alternative to improve the description of breakup is the use of hybrid 
models, in which hydrodynamic evolution in the early stage is combined with 
cascading in the late stage. In Figure 16, we compare results of one such hybrid 
model, Ultrarelativistic Quantum Molecular Dynamics (URQMD) (109; 110), to 
results of hydrodynamic calculations. Similar results were obtained by Teaney 
et al. (Ill; 112). Compared to hydrodynamic simulations, hybrid descriptions 
do not seem to notably reduce the overpredicted lifetimes. They tend to emit 
most of the pions at a time near or above 15 fm/c and significantly overpredict 
-Rout/-Rside ratios, whereas blast-wave parameterizations favor breakup times near 
or slightly below 10 fm/c. However, firm conclusions that the relative failure 
of hybrid models derives from the chosen equations of state cannot be made 
until comparisons are made between Boltzmann and hybrid calculations that use 
the same equation of state. Until such an analysis is performed, other issues 
will cloud the interpretation, such as whether viscous effects or details of the 
hydrodynamic/Boltzmann interface dominate results and might even invalidate 
the hydrodynamic approach. 

Entropy and pressure are intimately related, in that knowing the entropy den- 
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sity as a function of the energy density determines the pressure as a function of 
energy density. Whereas one can intuitively understand why a lower pressure 
would lead to longer lifetimes and larger source dimensions, the manifestations of 
changing the entropy are less transparent. One way to understand the effects of 
changing the entropy is to associate higher entropy with a larger effective number 
of light degrees of freedom. For instance, if one were to hadronize via a one-to-one 
parton-to-pion scheme, the volume per identical particle would change by a factor 
of the number of degrees of freedom. If the r.m.s. momentum is the same before 
and after hadronization, this would imply a change in entropy per particle equal 
to the logarithm of the ratio of effective degrees of freedom, and if ~ 40 light 
partonic degrees of freedom were immediately replaced by three pionic degrees of 
freedom, the system would lose ln(40/3) units of entropy per particle. To conserve 
entropy, a system must expand its volume by a similar ratio which implies an 
increase in radius parameters. Hydrodynamic models, which manifestly conserve 
entropy, use the energy stored in the latent heat to provide the heat necessary to 
preserve entropy during hadronization. Entropy conservation is more difficult to 
enforce in microscopic approaches. This underlines the challenges involved in ap- 
plying microscopic simulations in an environment of strongly interacting matter 
with ill-defined degrees of freedom, and it re-emphasizes the importance of gain- 
ing a better understanding of the validity of hydrodynamic calculations in such 
manifestly finite situations. It is peculiar that the entropy extracted from source- 
size and spectral measurements is consistent with the lattice-inspired equation 
of state (111; 112), whereas the source sizes extracted from hybrid models incor- 
porating similar equations of state significantly overpredict source sizes. Part of 
this contradiction can be explained by increases in the populations of baryons, 
which, because they inherently have more entropy per particle than do pions, can 
account for much of the missing entropy (114). Thus, the HBT puzzle does not 
necessarily imply an entropy puzzle. 

The HBT puzzle is not so much that radius parameters cannot be fit by models, 
but that our most sophisticated models, which incorporate a phase transition, fail 
to reproduce the data. The very gradual evolution of extracted source sizes as 
beam energies traverse a large range of energies is remarkable and puzzling in its 
own right. A simple explanation is that the equations of state do not dramatically 
change as the energy density changes from hadronic to super-hadronic densities, 
i.e., there is not even a hint of a latent heat. In fact, the failure of hydro dynamical 
models to reproduce femtoscopic measurements might be largely due to unrealis- 
tically large latent heats needed to reproduce elliptic flow data (115). However, 
as emphasized above, a host of unresolved issues prevent more quantitative con- 
clusions. These qualifiers cannot be lifted until much more thorough analyses of 
models are performed that entail a systematic exploration of the sensitivity of 
model predictions to both parameters and assumptions. This would necessitate a 
tremendous commitment from the community, but without it, many conclusions 
about the matter created in relativistic heavy ion collisions will remain vague. 

6 Summary 

Twenty-five years ago, the goal of femtoscopy was to demonstrate that one could 
measure a hadronic length scale with correlations, and if a result was on the order 
of a few Fermi, the analysis was deemed a success. In contrast, with the improved 
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accuracy of measurements, the enormous increase in statistics, and the simulta- 
neous development of phenomenology and theory, femtoscopy is now considered a 
precision measurement. Ten percent deviations between theory and experiment 
are now taken seriously as evidence that the spatio-temporal description of a 
model is significantly flawed. At RHIC energies, all six dimensions of the correla- 
tion function have been exploited to provide truly three-dimensional insight into 
the phase space cloud for particles of any momentum with any direction. 

For relativistic heavy ion collisions, there were expectations that a transition 
from hadronic to partonic matter might be accompanied by a large latent heat, 
which would bring about a dramatic change in the dynamics as beam energies 
traversed the range for exploring the mixed phase. Within this range, it was 
expected that the latent heat and the associated softening of the equation of 
state would manifest itself by slowing the explosion with lifetimes approaching 
or exceeding 20 fm/c. The signal of the phase transition would have been an 
increase of the effective lifetime for a range of beam energies followed by a return 
to more explosive and shorter-lived reactions at even higher energies. 

Extended lifetimes were not observed. Increasing beam energies from AGS to 
RHIC indeed causes larger energy densities and higher multiplicities, which push 
toward increasing the source volumes. However, much of this increase in multi- 
plicities is absorbed by higher phase space densities rather than larger source sizes. 
Combined with the increasing strength of radial flow, which provides smaller re- 
gions of homogeneity relative to the overall source volume, the result is that the 
effective dimensions change remarkably little over a wide range of beam energies. 
Furthermore, it appears that lifetimes of the reaction never leave the neighbor- 
hood of 10 fm/c. Not only does this represent a lack of evidence for a latent heat, 
it also suggests that there is no such latent heat. 

These conclusions remain only modestly guarded. Theory and phenomenol- 
ogy are progressing, but improvements in such aspects as mean-field effects or 
accounting for the smoothness approximation are not expected to change con- 
clusions by more than 10%. As with the improvements in including Coulomb 
effects in Section 4, removing some of the distortions and aberrations from the 
analyses is likely to be significant for fine-tuning models, but should not alter the 
conclusion that there is no large latent heat associated with the reaction. 

The theory of modeling, i.e., generating the source functions, has the greatest 
need for progress, because femtoscopic measurements are determined solely by 
the geometry of breakup, changes in chemical or kinetic evolution may have 
a significant impact on extracted source dimensions. The next generation of 
transport theories should be more flexible and will probably incorporate numerous 
effects such as in-medium mass changes, in-medium reduction of scattering cross- 
sections, viscous effects, and dynamical solutions for chemical rates. 

Despite the progress listed above, measurement has only begun to address 
the rich expanse of information available in correlations. Nearly all the three- 
dimensional analyses have been focused on identical-pion correlations. The huge 
data sets of the recent and upcoming runs at RHIC make it possible to analyze 
source functions for many pairs of particles in six full dimensions. In addition to 
providing important verification of identical-pion measurements, these analyses 
address other issues, such as whether all species flow and break up together. 

This is not a field for the complacent. As emphasized above, efforts at RHIC 
are just beginning to explore wholly new classes of correlations. Energy densities 
at the LHC might surpass those at RHIC by the same factor that those at RHIC 
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surpassed the AGS. Just as our visions of the future from twenty years ago proved 
largely naive, we should be prepared to be surprised with the femtoscopy of the 
next 25 years. 
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